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Abstract 

This  thesis  considers  the  design  of  a  linear  quadratic  Gaussian  (LQG)  controller 
for  a  ground-based  adaptive-optics  telescope.  The  incoming  aberrated  image  is  reflected 
from  a  97-actuator  deformable  piezoelectric  mirror,  then  measured  with  a  Hartmann- 
type  wavefront  sensor.  A  Kalman  filler  processes  the  outputs  of  the  wavefront  sensor  and 
obtains  estimates  of  system  states.  A  linear-quadratic  (HQ)  regulator  processes  these  state 
estimates  and  determines  an  appropriate  set  of  commands  for  the  deformable  mirror. 

Atmospheric  distortion  is  modeled  as  a  set  of  fourieen  Zcrnike  coefficients  whose  dy¬ 
namic  behavior  is  produced  by  excitation  of  a  set  of  shaping  filters  by  zero-mean  Gaussian 
white  noise.  The  response  of  the  mirror  to  control  voltages  is  modeled  as  a  set  of  Zernike 
coefficients  whose  dynamics  are  modeled  as  deterministic  first-order  systems.  The  entire 
control  system  is  simulated  using  the  Multimode  Simulation  for  Optima!  Filter  Evaluation 
(MSOFE)  software. 
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DESIGN  OF  A  LINEAR  QUADRATIC  GAUSSIAN  CONTROL  LAW  FOR  AN 

ADAPTIVE  OPTICS  SYSTEM 


1.  Introduction 


1.1  Motivation 

Certain  aspects  of  the  Air  Force  mission  require  high-resolution  imaging  of  distant 
objects.  The  resolution  attainable  is  theoretically  diffraction  limited  by  the  receiving  aper¬ 
ture.  From  physics,  this  diffraction  limited  resolution  is  a  function  of  both  the  diameter 
of  the  aperture  and  the  wavelength  of  light.  The  larger  th„  ratio  of  aperture  diameter  to 
wavelength,  the  better  the  resolution.  When  some  or  all  of  the  space  between  the  object 
and  aperture  contains  turbulent  atmosphere,  as  with  ground-based  imaging  of  celestial 
bodies,  the  resolution  is  degraded  considerably— the  image  wanders,  becomes  fuzzy  and 
distorted,  and  undergoes  intensity  fluctuations  [2].  It  has  been  reported  [l 3:360]  that,  on 
the  average,  these  effects  may  degrade  resolution  of  larger  telescopes  to  two  arc  seconds 
or  more,  corresponding  to  diffraction-limited  viewing  through  only  a  6-cm  aperture.  One 
obvious  approach  to  overcoming  atmospheric  distortion  is  to  place  the  telescope  above  the 
atmosphere.  An  example  of  this  approach  is  the  Hubble  space  telescope. 

The  Air  Force  Weapons  Laboratory  (AFWL)  at  Kirtland  AFB,  NM,  the  Rome  Air 
Development  Center  (RADC)  at  Griffis  AFB,  NY,  and  a  number  of  other  Government 
and  private  institutions  are  interested  in  compensating  for  these  deleterious  effects  by 
using  adaptive  optics  in  the  receiving  optical  system.  Adaptive  optics  is  the  use  of  active 
optical  components  for  compensation  of  unwanted  time-varying  optical  characteristics. 
For  example,  an  automatic  focus  mechanism  in  a  35-nun  camera  is  an  adaptive  optics 
system.  In  particular,  the  AFWL  is  actively  pursuing  the  use  of  adaptive  optics  using 
deformable  mirrors  for  ground-based  observation  of  orbiting  satellites,  and  is  the  sponsor 
of  this  research. 
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1.2  Research  Objective 

Tlie  overall  objective  of  this  research  is  to  develop  a  Linear  Quadratic  Gaussian 
(LQG)  control  law  for  an  adaptive  optics  system.  Such  a  control  law  provides  optimal 
control  in  the  sense  of  minimizing  a  defined  cost  function  associated  with  deviation  from 
desired  behavior.  LQG  also  accounts  for  uncertainty  in  measurements  and  the  stochastic 
nature  of  the  state  dynamics.  A  detailed  discussion  of  LQG  control  is  available  in  the 
literature  [25).  To  achieve  the  LQG  design,  several  tasks  are  required: 

1.  Define  a  state  space  on  which  to  base  design  and  analysis. 

2.  Develop  the  stochastic,  equations  which  model  the  open-loop  system. 

3.  Develop  the  LQG  control  law. 

4.  Simulate  and  evaluate  the  controlled,  closed-loop  system. 

The  first  task  involves  determining  a  method  of  quantifying  distortion  of  a  time- varying  im¬ 
age.  The  method  should  account  for  most  of  the  distortion  using  as  few  states  as  possible. 
The  states  should  also  be  conducive  to  measurement.  Once  the  states  are  identified,  the 
second  task  is  to  represent  the  stochastic  nature  of  their  dynamics  and  measurements  in  a 
&et  of  equations  (models)  which  describe  the  random  nature  of  the  processes.  Key  compo 
nents  of  such  a  model  arc  atmospheric  distortion,  mirror  response,  and  wavefront  sensing. 
Wherever  practical,  the  research  uses  actual  data  from  the  deformable  mirror  apparatus 
currently  in  the  optics  development  laboratory  at  the  AFWL.  Once  these  open-loop  “truth” 
models  are  developed,  the  actual  controller  can  be  designed.  This  task  includes  definition 
of  an  appropriate  cost  function,  design  of  a  Kalman  filter — possibly  of  reduced  order — to 
estimate  system  states  and  covariances,  and  design  of  a  linear-quadratic  (LQ)  regulator. 
The  final  task  includes  simulation  of  the  complete  system.  This  will  involve  implementing 
the  truth  models,  Kalman  filter,  and  controller  in  sofware  capable  of  simulating  random 
processes.  Multimode  Simulation  for  Optimal  Filter  Evaluation  (MSOFE)  [2Sl  is  such  a 
package. 
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1.3  Assumptions 

This  research  is  a  proof-uf -concept.  Several  assumptions  are  necessary  in  order  to 
carry  it  out  in  the  slotted  time.  Sonm  of  the  assumptions  are  made  based  on  literature 
precedent.  Many  of  the  assumptions  are  made,  to  account  for  lack  of  data.  Other  assump¬ 
tions  are  engineering  judgement. 

The  first  assumption  is  that  the  telescope  is  ground  based,  and  that  the  application 
is  satellite-viewing.  A  second  assumption  is  that  a  separate  control  system  is  responsible 
for  keeping  the  teloscop-  pointed  at  the.  sa.tellitt.  It  is  assumed  the  telescope  is  using  a 
visible  point  source  of  light  at  or  near  the  satellite  as  a  reference,  and  that  anisoplanatism 
is  negligible,  i.e.,  the  light  from  the  reference  traverses  the  same  atmosphere  as  the  light 
from  the  satellite.  This  implies  that  deforming  the  mirror  to  improve  the  reference  image 
will  also  sharpen  the  satellite  image.  The  wavelength  A  of  the  reference  image  light  is 
assumed  to  be  514  nanometers  (nm),  the  value  used  during  laboratory  evaluation  of  the 
hardware  [22]. 

The  deformable  mirror  is  a  97-clement  piezoelectric  mirror  whose  actuators  are  ar¬ 
ranged  in  evenly-spaced  rows  and  columns.  Furthermore,  based  on  limited  data,  it  is 
assumed  the  actuators  are  linear,  and  that  all  actuator  influence  functions  are  identical 
and  symmetric.  The  effects  of  additional  factors  (lenses,  vibrations,  etc.)  are  not  modeled. 
In  an  actual  telescope,  separate  active  tilt  mirrors  may  precede  the  deformable  mirror. 
These  non-deformable  tilt  mirrors  remove  the  gross  tilts  of  the  overall  image.  This  re¬ 
search  does  not  explicitly  model  these  mirrors,  but  assumes  they  remove  95  %  of  the  gross 
tilts. 

Although  the  statistics  of  atmospheric  distortion  of  images  are  time  varying — for 
example,  there  is  generally  less  distortion  at  night  than  day — the  stochastic  model  for  the 
atmosphere  is  assumed  time-invariant.  This  assumption  is  valid  if  the  duration  of  the 
observation  is  short  relative  to  the  time- varying  nature  of  the  atmospheric  statistics. 

Finally,  the  very  nonrestrictive  assumption  is  made  that  the  time  del.ay  due  to  light 
propagation  through  the  optical  components  is  negligible.  The  time  delay  due  to  finite 
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sampling  time,  however,  is  not  negligible.  The  effect  of  time  delay  in  the  sampling  process 
will  not  be  specifically  modeled,  but  will  be  used  as  an  upper  limit  on  the  sampling  rate. 

l.J  Treatment 

This  thesis  presents  material  in  a  topical  fashion.  First,  Chapter  2  begins  by  laying 
the  groundwork  fundamental  to  the  understanding  of  atmospheric  turbulence  and  its  effect, 
on  optical  viewing.  A  qualitative  discussion  of  an  adaptive  optics  system  is  presented.  The 
remaining  portion  of  Chapter  2  defines  the  state  space  used  for  design  and  analysis  of 
the  control  system.  Chapter  3  develops  stochastic  models  ol  the  atmospheric  distortion, 
deformable  mirror,  and  wavefront  sensor.  Chapter  4  presents  the  design  of  the  LQG  control 
law,  and  Chapter  5  discusses  the  closed-loop  simulation  and  evaluation  thereof.  Chapter 
6  presents  conclusions  of  this  research  and  recommendations  for  future  efforts. 


//.  Background 

2.1  \tmosjdn  "ic  Optic* 

.  .1.1  Introdurtum  The  atmosphere  can  have  several  effects  on  light  traversing  it. 
Certain  chemical  <  mstit  aents  of  the  atmosphere  may  selectively  absorb  certain  lreqnoncies 
of  light  am!  retransmit  the  energy  at.  different  frequencies.  If  the  light  is  a  high-energy  beam 
such  as  a  high-energy  la- er,  it  is  possible  for  the  light  energy  absorbed  by  tin  atmosphere 
to  affect  the  index  of  refraction  in  the  immediate  vicinity  of  the  beam,  v.j.ich  in  turn 
affects  the  shape  of  the  beam  profile.  This  phenomenon  is  known  as  blooming  [41:223]. 
In  the  presence  of  aeroso.s,  light  may  experience  spatial  discontinuities  in  the  direction 
of  propagation,  a  phenomenon  known  as  scattering.  Finally,  in  imaging  applications,  the 
turbulence  of  the  atmosplmro  can  cause  continuous  spatial  and  temporal  variations  in  the 
index  of  refraction  along  the  ray  path,  resulting  in  temporal  and  spatial  modulation  of 
the  received  intensity  and  phase.  The  modulation  of  intensity  is  known  as  scintillation 
[5:37],  an  example  of  which  is  the  twinkling  of  stars  [37:224]  [36:3].  As  a  rule  of  thumb, 
scintillation  effects  can  be  considered  insignificant  when  the  wavelength  A,  propagation 
distance  through  turbulence  /,,  and  telescope  entrano  aperture  diameter  D.  are  related 
by  [44:2580]  [45:819]  : 

y/IZaD  (2.1) 

For  example,  for  an  atmospheric  propagation  distance  of  20,500  meters  ]3S:A4],  wavelength 
of  514  nanometers  and  a  telescope  diameter  of  1  meter,  the  result  is  a  less  severe  inequality 
0.102  <  1.0.  Nevertheless,  much  of  the  literature  [9:1435]  [29]  [47:3]  [1C]  suggests  the 
intensity  modulation  is  negligible.  The  modulation  of  phase  is  the  leading  contributor  to 
image  wandering,  fuzziness,  and  distortion.  This  research  is  therefore  limited  to  correction 
of  phase. 

2.1.2  Turbulence  The  physical  property  cf  an  air  p;  .cel  which  most  influences 
phase  of  propagating  light  is  its  index  of  refraction.  Tins,  in  turn,  is  a  function  of  various 
physical  parameters  such  as  temperature,  pressure,  humidity,  ami  wavelength,  to  name  a 
few.  Empirical  equations  describing  the  dependence  of  the  index  of  refraction  on  these 
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parameters  are  generally  available  [5:10]  [18:531]  [31:101].  In  the  turbulent  atmosphere 
these  physical  parameters  are  generally  functions  of  time  and  space,  and  hence  index  of 
refraction  is  also  a  function  of  time  and  space.  Because  of  the  random  nature  of  turbulence, 
there  is  no  deterministic  expression  of  this  relationship. 

A  popular  representation  [5:12]  [18:330]  of  atmospheric  turbulence  considers  the  at¬ 
mosphere  to  consist  of  many  parcels  or  “eddies”  of  air  moving  with  respect  to  each  other. 
The  sizes  of  these  eddies  change  with  time;  the  larger  ones  generally  dissipate  into  smaller 
ones.  Eventually,  the  eddies  become  small  enough  such  that  their  kinetic  energy  dissipates 
into  heat.  The  sizes  of  the  eddies  are  characterized  relative  to  Lq  and  1 0,  quantities  having 
units  of  length  and  known  as  the  “outer  scale”  and  “inner  scale”  of  turbulence,  respec¬ 
tively.  Eddies  larger  than  the  outer  scale  are  anisotropic  and  are  formed  by  wind  shear 
and  temperature  gradients.  An  anisotropic  eddy  has  spatial  statistics  which  depend  on 
direction  [18:338].  These  eddies  are  said  to  lie  in  the  “input  range”  since  they  represent 
newly-formed  eddies,  inputs  to  atmospheric  turbulence.  For  eddies  smaller  than  Lo  but 
larger  than  1Q,  kinetic  energy  effects  are  more  significant  than  viscosity  effects,  and  the 
turbulence  is  essentially  isotropic.  An  isotropic  eddy  has  spatial  statistics  which  do  not 
depend  on  direction.  Eddies  of  this  size  are  sajd  to  be  in  the  “inertial  subrange”  since 
oertia  is  significant  relative  to  friction.  When  the  size  of  the  eddy  is  less  than  I0,  viscosity 
effects  dominate  over  kinetic  effects,  and  the  eddy  is  said  to  be  in  the  “dissipation  range” 
[18:336] .  At  low  altitude  Lo  is  on  the  order  of  meters,  whereas  Iq  is  on  the  order  of  millime¬ 
ters.  Both  the  local  inner  and  outer  scales  of  turbulence  generally  increase  with  altitude 
[5:12].  Light  propagating  vertically  through  the  atmosphere,  therefore,  encounters  eddies 
of  various  sizes  and  of  various  subranges.  The  end  result  is  a  randomly-distorted  image 
whose  degree  of  distortion  is  not  known  without  measurement. 

2.1.3  Structure  Functions  The  statistics  of  atmospheric  turbulence  are  often  rep¬ 
resented  as  “structure  functions”  [37],  The  presentation  of  the  structure  function  here 
closely  follows  that  of  Tatarski  [37]  and  Ishimaru  [18].  Any  direct  quotes  are  from  tire 
latter  source. 


Some  physical  characteristics  associated  with  a  given  location  in  the  turbulent  atmo¬ 
sphere  can  be  described  as  random  processes.  For  example,  let  h(t)  be  a  random  process 
describing  the  absolute  humidity  at  a  specified  location.  This  random  function  . .  is  not 
strictly  stationary.”  The  statistics  of  h{()  can  and  generally  do  change  with  time.  The 
difference  h(t  +  r)  —  h(i),  however,  is  considered  stationary.  This  property  is  accurate 
for  many  atmotphcric  variables.  Thus  it  can  be  said  the  function  h(t)  has  stationary 
increments: 

£(h(t  +  r)-/i(0W(r)  (2.2) 

where  t  [•]  represents  the  expectation  operator  and  /(r)  is  a  function  of  the  time  difference, 
t  .  The  temporal  structure  function  for  this  example,  D(r)  ,  is  defined  as: 

Vh(r)  =  £  [|M<  +  r)  -  /r(0|2]  (2-3) 

The  correlation  function  is  a  more  familiar  statistical  relation: 

Bn(tnt2)  =  £[h{h)h{t2))  (2-4) 

and  is  related  to  the  structure  function  by: 

Vh{T)  =  6h(t  +  r,t  +  r)  T  -  Bh(t  +  r,t)  -  Bh{t,t  4-  r)  (2.5) 

The  spatial  analog  to  the  random  process  with  stationary  increments  is  the  “locally 
homogeneous”  random  function.  For  example,  let  h{ r)  be  a  random  “process”  describing 
the  absolute  humidity  at  a  specified  time  at  location  r  .  For  convenience  the  /i(-)  function 
symbology  is  reused  here.  This  random  function  “. . .  is  not  strictly  homogeneous.”  The 
statistics  of  h( r)  can  and  generally  do  change  with  location.  The  difference  h(r+p  )  -  h( r), 
however,  is  homogeneous.  Thus  the  function  h( r)  is  said  to  be  locally  homogeneous.  The 
spatial  structure  function  for  this  example,  V(p)  ,  is: 

Vh{p)  =  S  |Mr  +  £)-/;(r)|2  f.2.6) 
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where  the  !)/,(■)  symbology  is  reused  for  convenience.  If  the  process  is  isotropic  as  well  as 
locally  homogeneous,  the  spatial  structure  function  depends  only  on  the  magnitude  of  tie- 
spatial  separation: 


p)  =  £ 


h{r  +  p)  -  A(r)j 


P  = 


(2.7) 


A  key  assumption  in  applying  the  temporal  and  spatial  structure  functions  to  atmo¬ 
spheric  turbulence  is  the  concept  of  “frozen  turbulence"'  [5:17],  also  known  as  “Taylor’s 
hypothesis”  [19:133],  This  assumption  implies  the  temporal  fluctuations  in  various  mete¬ 
orological  variables  at  a  given  location  are  caused  by  a  “snapshot”  realization  of  spatial 
fluctuations  flowing  by.  Taylor's  hypothesis  is  “usually  a  good  approximation  for  optical 
propagation  through  the  atmosphere.”  [19:133] 


.2.1.4  Kolmogorov  Statistics  The  index  of  refraction  n  at  a  location  r  can  be  ex¬ 
pressed  as  the  sum  of  the  mean  and  a  fluctuation  about  the  mean  [5:17]: 


n( r)  =  S  (n(r)]  +  ni(r) 


(2.8) 


where  n^r)  is  a  zero-mean  Gaussian  random  variable.  The  spatial  structure  function  for 
:h-  isotropic  random  variable  describing  the  index  of  refraction  fluctuation  is: 


vni{p)  =  £ 


ni(r+  p)  -  ni(r)! 


(2.9) 


The  unity  subscript  in  Equation  (2.9)  is  often  omitted,  with  the  Auction  definition  implied. 

Kolmogorov’s  famous  result,  a.s  described  by  several  references  [6:155]  [2:11]  [17:22] 
[31:101]  [34:288]  [5:19]  [19:526]  [9:1430],  is  that,  for  spatial  separation  within  th  >  inertial 
subrange  t'o  p  <C  To,  Equation  (2.9)  is  closely  approximated  by: 


Vn(p)  =  C-p 


-  r-n 2/3 


(2.10) 


where  is  the  refractive  index  structure  constant,  a  parameter  which  indicates  the  degree 
to  which  atmospheric,  turbulence  affects  optical  propagation  [34:2S8j.  C %  is  generally  mod¬ 
eled  as  a  function  of  altitude  and  time  of  day  [5:20];  several  models  describe  it  as  functions 


Figure  2.1.  Refractive  Index  Structure  Constant 

thereof  [13:290]  [5:20]  [6]  [37:A4]  [30:102].  Two  such  models  are  the  SLC-Day  model  and 
the  Hufriagel-Stanley  model  [30:102].  Figure  2.1  shows  that  these  widely-accepted  models 
are  in  less  than  perfect  agreement,  and  that  the  altitude  dependency  of  C ^  varies  over  three 
orders  of  magnitude  [30] . 

Since  the  phase  of  a  wavefront  is  of  interest  for  many  adaptive  optics  applications, 
one  may  wonder  if  it  has  a  structure  function.  Fried  [9:1430]  discusses  the  phase  structure 
function: 

'£*(/»)  =  £  (jp(l  +  £)  -<?(r)|2  (2.11) 

which  is  actually  defned  as  a  phase  fluctuation  structure  function,  since  p  is  the  deviation 
of  phase  about  the  aperture-averaged  phase.  Based  on  the  Kolmogorov  turbulence  model 
and  neglecting  intensity  fluctuations,  Fried  writes: 

iyvp)  =  Aps'3  (2.12) 


n.  ^ 
i.'  o 


where  the  parameter  A  depends  on  propagation  path,  wavelength,  and  environmental 
conditions.  He  then  defines  the  coherence  length  parameter  ro  in  terms  of 'A: 

r0  ^  (G.88/^)3/5  (2.13) 

and  thus: 

V^p)  =  G.88  (p/r0)5/3  (2.14) 

This  ro  parameter  is  generally  known  as  Fried’s  coherence  length  in  the  later  literature 
[29:209]  [11:550]  [40:1774]  [17:159G]  [31:102].  He  describes  it  as  the  “diameter  of  a  hetero¬ 
dyne  collector  for  which  distortion  effects  begin  to  seriously  limit  performance”.  Another 
interpretation  is  that  it  is  the  diameter  of  an  aperture  such  that  the  rms  phase  distortion 
is  1  radian  [23].  He  finishes  his  introduction  by  stating  that,  for  visible  and  near-lit  prop¬ 
agation  from  an  approximately  zenith  source,  typical  values  of  r0  are  on  the  order  of  a  few 
centimeters. 

An  equation  for  the  coherence  length  in  terms  of  other  parameters  is  given  by  Parenti 
[31:102]: 

/  a2  X3/5 

r°  =  °'185  (iic(C)  JdhCl(h))  (2,13) 

where 


To 

Fried’s  coherence  length  (in) 

A 

wavelength  (m) 

c 

zenith  angle  of  source  (rad) 

clih)  = 

refractive  index  structure  constant  (m-2/3) 

h  '  = 

altitude  (nr) 

Wang  and  Markey  [46:78]  give  a  simplified  equation  for  r0  where  propagation  is  from 
the  zenith  and  the  refractive  index  structure  constant  is  assumed  to  be  constant  along  the 
path: 

r0  =  1.68(C2  z  fc2)-3/5  (2. 1C) 

where 
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z  —  path  length  through  turbulence  (m) 

k  —  wavenumber  =  2 /" / A  (m_1) 

£  =  ?.cnith  angle  of  source  (rad) 

Walters  and  others  [42:8281  present  some  measured  values  of  ro  for  vertical  paths  at 
mountaintop  sites  in  the  White  Sands,  New  Mexico  region.  The  nighttime  average  (over  a 
seven-month  period)  was  found  to  be  'J.O  cm  ±  4.0  cm  (1-ct);  the  daytime  average  was  4.5 
cm  ±  1.8  cm  (l-o).  Starlight  was  used  as  the  light  source. 

Parent';  [31:104]  gives  some  estimates  of  r0  derived  from  the  liufnagel-Valley  and 
SLC-Day  models  of  C^(h)  at  500-nm  wavelength  and  two  zenith  angles.  These  are  shown 
in  Table  2.1. 


Table  2.1.  Calculated  Coherence  Length 


Model 

Zenith  Angle 

t-o 

1IV-21 

0° 

5.0  cm 

HV-21 

45° 

4.0  cm 

SLC-Day 

0° 

5.1  cm 

SLC-Day 

45° 

4.1  cm 

The  significance  of  the  development  presented  thus  far  is  that  the  effect  of  atmo¬ 
spheric  turbulence  on  optical  viewing  is  complex.  The  time  behavior  of  the  quality  of 
the  received  image  depends  on  many  factors  such  as  the  wavelength,  altitude-dependent 
structure  function  constant,  zenith  angle,  and  altitude  of  the  telescope.  The  wind  velocity 
and  its  altitude  distribution  also  affect  image  quality'.  In  addition,  satellite  motion  with 
respect  to  the  observer  is  a  factor.  The  introductory  discussion  presented  thus  far  could 
continue  to  the  point  of  having  an  analytically-derived  statistical  model  of  the  effect  of 
atmospheric  turbulence  on  image  quality,  but  the  mathematical  rigor  is  beyond  the  scope 
of  this  research.  The  reader  is  referred  to  the  literature  for  an  appreciation  of  the  rigor 
involved  [16].  The  significant  result  is  that  no  single  model  will  be  accurate  for  all  possible 
atmospheric  conditions.  The  use  of  multiple  models  will  be  briefly  discussed  in  Chapter 
VI. 


2.2  Zemike  Functions 


2.2.1  Introduction  It  is  obviously  necessary  to  have  some  means  of  expressing  the 
phase  distortion  present  in  an  i  mage.  This  section  develops  the  use  of  the  Zernike  basis 
functions  as  such  a  means. 

Assume  a  circular  aperture  of  diameter  D  is  pointed  at  a  coherent,  monochromatic 
point  source  of  light  located  a  distance  L  away.  Further  assume  that  L  is  much  greater 
than  D.  If  there  are  no  distortive  elements  along  the  propagation  path,  the  phase  of  the 
incident  wave  within  the  aperture  would  be  spatially  constant,  i.e.,  the  phase  of  any  two 
points  within  the  aperture  would  be  the  same.  It  could  also  be  said  for  such  a  case,  that 
the  phase  at  any  point  within  the  aperture  is  equal  to  the  spatially  averaged  phase.  Figure 
2.2  shows  a  two-dimensional  representation  of  this  situation. 


If  a  turbulent  atmosphere  is  now  introduced  into  some  or  all  of  the  space  between 
source  and  observer,  the  wavelionts  at  the  aperture  will  be  distorted,  i.e.,  two  points  within 
the  aperture  will  not  necessarily  have  the  same  phase.  Furthermore  the  phase  at  any  point 
within  the  aperture  will  not  necessarily  equal  the  spatially  averaged  phase.  Figure  2.3 
shows  a  two-dimensional  representation  of  this  situation. 
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Figure  2.3.  Distorted  Wavefronts  at  an  Aperture 

The  absolute  phase  of  the  electromagnetic  wave  in  the  aperiuie  is  a  function  of  time 
and  of  the  spatial  (polar  or  rectangular)  coordinates  within  the  aperture: 


4>  =>  4>{r,Q,t) 
<p  =>  4>(  x,y,t) 


(2.17) 

(2.1C) 


where 


absolute  phase 
radial  distance 
angular  coordinate 
x-coordinate 
y- coordinate 


=  x1  +  jr 
=  tan_1(j//x) 

=  t  cos(0) 

=  r  sin(0) 


The  functions  represented  by  Equations  (2.17)  and  (2.18)  do  not  have  identical  functional 
forms;  both  functions  are  given  the  same  name  <t>  for  convenience.  Phase  can  be  expressed 
in  units  of  degrees,  radians,  or  wavelengths.  Likewise,  x,  y ,  and  r  can  be  expressed  in  units 
of  wavelengths  or  the  length  normalized  by  the  rail  us  of  the  aperture. 


The  right-hand  sides  of  Equations  (2.17)  and  (2.18)  are  expressible  in  several  ways. 
One  option  might  be  to  sample  the  phase  at  numerous  locations  in  the  plane  and  let  the 
vector  of  samples  define  the  shape  of  the  phase  function.  This  “vector  space"  approach 
has  been  used  by  AFIT  researchers  [3]  [12]  [27]  [33].  Another  approach,  the  one  taken  in 
this  research,  is  to  express  the  phase  as  a  sum  of  functions: 

N 

<Hr,&,t)=YlFi(r,Q,t)  (2.19) 


More  conveniently,  the  phase  can  be  expressed  as  a  linear  combination  of  basis  func¬ 


tions: 


N 


£ 


a,(0  Zi(r,0) 


(2.20) 


where  the  Z,(r,0)  are  basis  functions  which  span  the  functional  space  containing  d>(r,  0), 
and  the  a,(t)  arc  the  coefficients.  The  choice  of  which  set  of  basis  functions  to  use  is  a 
design  decision.  Some  mentioned  in  the  literature  include  Legendre  [35],  Karhunen-Loeve 
[29]  [43]  [46],  and  Zernike  [4]  [9]  [1G]  [46]  [44]  [7]  [45]  [15]  [38]  functions.  In  the  case  of 
turbulence  modeled  with  Kolmogorov  statistics,  the  Karhunen-Loeve  functions  are  not 
analytic  [29:210],  The  Zernike  set  of  basis  functions,  analytic  by  definition,  are  used  as 
basis  functions  in  this  research. 


2.2.2  Zernike  Functions  The  Zernike  functions  are  the  Z,(r,  0)  ol  Equation  (2.20). 
Each  Zernike  function  is  a  real-valued,  dimensionless,  deterministic  function  of  position 
within  the  aperture.  The  a,(t)  are  the  Zernike  coefficients.  Each  Zernike  coefficient  is  a 
real-valued  function  of  time,  with  units  of  wavelength  for  this  research.  The  i-th  Zernike 
function  can  be  expressed  as  the  product  of  a  radial  function  and  an  azimuthal  function 
[4]: 

Z,(r,©)= /,(r)y,(0)  (2.21) 

Details  of  generating  the  radial  and  azimuthal  functions  of  Equation  (2.21)  are  available 
in  the  literature  [38:79] [4].  Table  2.2  presents  the  first  fifteen  Zernike  functions  (0 — 14) 
along  with  their  radial  and  azimuthal  order,  n  and  in,  respectively.  The  functions  are  also 
expressed  in  rectangular  coordinates  as  well.  Three-dimensional  plots  of  the  first  fifteen 
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Table  '2.2.  Zornike  Functions 


m 

n 

in 

Zi[r.Q) 

0 

0 

0 

1 

1 

1 

1 

1 

2(£)cos(0) 

fa 

2 

1 

1 

2(fi)sin(0) 

fa 

■ 

2 

0 

^(2(-;-V*-l) 

#(2x’2  +  2i/  -  n 2) 

Bj 

2 

2 

v/G(^)‘sin(20) 

i ${xy) 

m 

2 

2 

\Aiifj)2  cos(2  Q) 

#(*’  -  r ) 

G 

3 

1 

V8(3(£)3-2(£))sin(0) 

${3x7  +  3y2  -  2 R2)y 

7 

3 

1 

x/8(3(^)3-2(^))cos(0) 

+  3?/*  —  2R2)i 

8 

D 

3 

\/$(^)3sin(30) 

-  y‘2)u 

9 

n 

3 

x/S(^)3cos(30) 

$(x2  -  3y2)j 

10 

4 

0 

^(C(^)4-C(^)2+l) 

)^(6(.t2  +  t/2)2  -  GR2(x2  +  t/2)  +  A'1) 

11 

m 

2 

V^(4(^)4-^(s)2)cos(20) 

12 

S 

2 

V^0(4(^)4  -  3(^)2)siit(20) 

2)7p(4z2  -•  4r/‘  31l2)xy 

13 

n 

4 

v/To(^)‘l  cos(40) 

^{x4  -  6 x2y2  +  y4) 

14 

n 

n 

x/lO(£)4sin(40) 

AfaP(x2  ir)xy 

Zcrnike  function  are  shown  in  Appendix  A. 

One  advantage  of  using  the  Zcrnike  functional  space  is  that  the  functions  correspond 
to  aberrations  commonly  studied  in  optics  [20:196].  Table  2.3  shows  the  names  of  the  more 
common  ones.  By  convention,  x-tilt  is  defined  in  this  research  to  be  tilt  about  the  x-axis. 

The  Zernike  functions  form  an  orthogonal  basis  set  which  satisfies  [38:79]: 

1  f2~  fR 

/  dQ  dr  r  Zi[r,0)  Zj(r,G)  =  tl}  (2.22) 
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Tabic  2.3.  Zcrnikc  Function  Names 


i 

Aberration 

0 

piston 

1 

y-tilt 

2 

x-tilt 

3 

focus 

4 

astigmatism 

5 

astigmatism 

0 

coma 

7 

coma 

where  <5,y  is  t.lie  Kronecker  delta: 


i 


i  =  ; 
*  #  i 


(2.23) 


As  shown  in  Appendix  C,  a  consequence  of  Equation  (2.22)  is  that  the  root-mean-square 
(rms)  value  of  phase  in  an  aperture  is  the  square  root  of  the  sum  of  squares  of  the  corre¬ 
sponding  Zernike  coefficients: 


~  \f  o-l  +  G1 


0-2  +  •  •  •  +  Ok 


N 


(2.24) 


The  expansion  coefficients  for  an  arbitrary  phase  associated  with  Equation  (2.20)  can  be 
obtained  by  [38:79]: 


.m  =  f  dQ  f  <*rr  lV(r,QMr,0,Ofr(r,0) 

1  f  dQ  f  dr  r  W(r,  0) 


where  W(r,0)  is  the  aperture  weighting  function  defined  by  [38:79]: 


W(r,0) 


✓ 

I  1  r  <  R 
jo  r  >  R 


(2.25) 


(2.26) 


2-12 


It  should  be  noted  that  this  definition  of  W(r,0)  is  selected  to  be  consistent  with  the 
scaling  coefficients  in  Z,(r,  0).  For  a  circular  aperture  of  known  radius,  the  phase  within 
the  aperture  at  a  given  time  can  also  be  expressed  in  shorthand  notation  as  a  vector  a  of 
expansion  coefficients: 


n1  ( t )  = 


«o(0  '-,(/)  a;(<) 


«A’(0 


(2.27) 


rn 

Post-multiplication  of  n  1  (/)  by  a  column  vector  of  A'  +  1  Zcrnike  functions  is  implied  in  this 
notation.  In  general,  N  =  oo  is  required  to  obtain  an  exact  description  of  an  arbitrary  phase 
front  in  terms  of  basis  functions.  Based  on  spatial  Nyquist  consideration  of  the  actuator 
and  wavefront  sensor  geometry  for  the  optics  system  at  the  AFWL,  and  of  relative  mean 
square  phase  error  content  of  each  Zcrnike  mode  in  the  presence  of  Kolmogorov-modeled 
turbulence  (see  Appendix  B),  the  truth  model  in  this  research  is  selected  to  consist  of  the 
first  fifteen  Zernike  modes  (K  =  14). 

As  an  example,  of  what  an  arbitrary  sum  of  the  first  15  Zernike  functions  would  look 
like,  Figure  2.4  shows  a  plot  of  Equation  (2.20)  for  a  coefficient  vector  of: 


rp  r 

ar(0  = 


20  1  1  1  1  1  1  3  1  1 


1  1 


(2.28) 


The  zeroth  Zernike  mode  (i.e.  piston)  is  not  a  distortive  contributor.  Also,  it  is  not 
measurable  by  the  sensor  used  to  measure  wavefront  distortion.  Thus  we  can  subtract  the 
piston  contribution  from  the  absolute  phase  and  obtain  a  functional  space  which  describes 
the  deviation  in  phase  from  the  aperture-average  phase.  From  this  point  forward,  only 
Zernike  modes  1-14  are  considered,  and  <^»(r,  0,  t)  is  redefined  to  be  the  deviation  from 
aperture-averaged  phase,  also  termed  phase  distortion.  Thus,  using  the  shorthand  notation 
introduced  earlier: 


<j>{r,0,t)  =>  ar(f)  = 


fl](0  oz(0  ••• 


(2.29) 


Figure  2.5  shows  the  same  plot  as  Figure  2.4,  but  with  the  piston  (i.e.  average  phase)  not 
included. 
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Figure  2.4.  Weighted  Sum  of  Zernike  Functions  0-14 

One  final  topic  is  required  to  be  discussed  regarding  Zemikc-function  description  of 
optical  phase  deviation — scaling.  The  question  is,  if  a  large  aperture  sample  of  an  incoming 
wavefront  is  focused  into  a  smaller  aperture,  what  is  the  effect  on  the  expansion  coefficients? 
The  answer  is:  nothing.  For  example,  consider  the  geometry  of  Figure  2.6.  Aberrated  light 
from  a  distant  source  enters  a  ono-meter  radius  circular  aperture.  Assume  the  aberration 
consists  of  ten  wavelengths  of  x-tilt  across  the  aperture.  At  this  large  aperture: 

2 

<t>{x,y)  =  «2  Z2{x,y)  =  a2  —  y  =  2  a2  y  (2.30) 

At  y  =  1  the  phase  deviation  is  five  wavelengths;  therefore  a2  —  5/2  for  the  large 
aperture. 

By  passing  through  a  set  of  perfect  lenses,  the  size  of  the  image  is  reduced  lo  a  circle 
of  0.1 -meter  radius.  At  this  small  aperture: 

0\z>y)  ~  a2  7\ [x,y)  -  a2  ~y  =  20  a'2  y  (2.31) 
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figure  2.o.  Weighted  Sum  ofZeinike  Functions  1-14 

At  y=0.1,  the  phase  deviation  is  still  five  wavelengths;  therefore  a'2  =  5/2  for  the  small 
aperture  as  well.  This  same  reasoning  can  be  extended  to  all  the  Zernike  modes.  This 
aperture  invariance  is  convenient  for  imaging  applications,  since  generally  a  large- aperture 
telescope  gathers  light,  a  medium-aperture  deformable  mirror  processes  it,  and  a  small- 
aperture  wavefront  sensor  measures  its  distortion. 

2.3  Stochastic  Slate-Space  Modeling 

The  development  of  a  Linear  Quadratic  Gaussian  (LQG)  control  law  requires  a  linear 
stochastic  state-space  model  of  the  system  of  interest.  The  required  forms  for  continuous 
dynamics  and  discrete  measurements  aie  [24:169]: 

s(0 

£(*.) 

where 


=  F(t)x(l)  +  B(t)u(r)  +  G(t)w{t)  (2.32) 

=  H  ltt)x[U)  +  ylti)  (2.33) 
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Figure  2.6.  Change  in  Radius 

time  rate-of-change  of  system  state  vector 

system  state  vecioi 

state  dynamics  matrix 

deterministic  input  vector 

input  distribution  matrix 

zero-mean  white  Gaussian  driving  noise  vector 

driving  noise  distribution  matrix 

discrete  measurement  vector 

discrete  output  matrix 

zero-mean,  discrete-time  white  Gaussian  measurement  noise 


The  strength  of  the  zero-mean  white  Gaussian  driving  noise  w(i)  is  defined  by 
[24:155]: 

fU(02£T(0}=Q(0«U-0  (2-34) 

where  Q (t)  is  the  noise  strength  matrix  and  6(  )  is  the  dirac  delta. 

The  covariance  of  the  zero-mean  measurement  noise  v(<,)  is  defined  by  [24:174): 

£{v(OvT(f;)}  =  R  (u)6,}  (2.35) 


where  R(t,)  is  the  noise  covariance  matrix. 
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The  dimensionalities  of  the  vectors  and  matrices  associated  with  continuous  dynamics 
and  discrete  measurements  are  given  in  Table  2 .4  [25:9]. 

Table  2.4.  Dimensionalities  of  State-Space  Matrices  and  Vectors 


Matrix/ Vector 

Dimension 

F«) 

7i  X  n 

x(0 

71  x  1 

B(0 

71  X  r 

u(0 

r  x  1 

0(0 

71  X  £ 

s(0 

JXl 

Q(0 

s  x  s 

z(b) 

m  x  1 

H(b) 

m  x  n 

x(b) 

m  x  1 

R(b) 

TO  X  777. 

If  the  dynamic  system  and  measurement  device  are  time-invariant,  which  was  as¬ 
sumed  for  the  nominal  models  of  this  research,  the  time-dependence  can  be  dropped  from 
the  notation  associated  with  the  defining  matrices: 


£(*) 

=  Fx(t)  +  Bu({)  +  Gw(() 

(2.36) 

z(b) 

=  Hx(b)  +  v(b) 

(2.37) 

The  dynamic  driving  noise  strength  and  measurement  noise  variance  lose  their  time- 
dependency  and  become  Qand  R,  respectively.  Such  time  invariance  implies  that  the 
statistical  behavior  of  the  system  does  not  change  with  time. 

Kalman  filters,  which  will  be  discussed  shortly,  often  model  the  continuous- time 
system  dynamics  of  Equation  (2.36)  in  the  discrete-time  domain.  An  '‘equivalent”  discrete- 


time  model  of  the  linear  system  dynamics  can  be  written  as  [24:171]: 


*(<«+ 1)  =  4>{ti+ 1  -  b)x(t,)+  Brfu(J,)  +  w^Q) 


(2.38) 


where  <f>(t,+1  -  t.)  is  the  n  x  n,  time-invaria.nl  state  transition  matrix  associated  with 
F  [24:41]  and  the  zero-mean  discretized  driving  noise  wj(/,)  has  strength  Q defined  by 
[24:171]: 


=  QAj  =  [/  + 
IA 


4>(t,+1  -r)GQGTchT(t1+1 


-  t)  dr 


(2.39) 


The  discrete-time  input  distribution  matrix  is  defined  by  [24:171]: 


(2.40) 


Chapter  III  presents  the  definitions  of  these  vectors  and  matrices  in  terms  of  the 
adaptive  optics  system  of  this  research. 


Linear  Quadratic  Regulator 

The  simple  LQG  regulator  concept  is  based  on  the  assumption  that  there  are  costs 
associated  with  nonzero  states  (x)  and  nonzero  control  inputs  (u).  Furthermore,  the  costs 
are  assumed  to  be  quadratic  in  nature,  that  is,  the  cost  is  proportional  to  the  weighted 
squares  of  the  states  and/or  control  inputs.  Maybeck  [25:10]  defines  a  simple  form  of 
quadratic  cost  function  as: 


J  =  £ 


[xT(t,)X(C)x(t,)  +  uT(/.,)U(t,)u(t,-)]  +  ^xT(!A'+1)Xfx(tA/+i) 

-_n  ^  ^ 


(2.41) 


where  to  and  <A+1  are  the  start  and  final  times,  respectively.  The  i  and  N  are  time  indices 
here  and  should  not  be  confused  with  their  use  as  mode  indices  in  Zernike  functions.  A 
more  general  cost  function  includes  cross  terms  between  states  and  control  inputs  [25:73j, 
but  this  extension  is  not  used  here. 

The  design  of  the  LQG  regulator  involves  determining  the  optimal  input  voltages  to 
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the  mirror,  u*(f,),  which  minimizes  the  right-hand  side  ofEquation  (2.41).  Minimization 
of  this  cost  must  be  related  to  minimization  of  the  phase  distortion  in  the  image  reflected 
from  the  deformable  mirror  (maximizing  performance).  This  relationship  is  implemented 
via  the  weighting  matrices  X(l,)  ,  U(t,)  ,  and  X  j  . 

Since  only  noise-corrupted  measurements  are  available,  a  Kalman  f'ter  is  used  to 
estimate  the  states  (x,)  given  the  measurements.  Figure  2.7  shows  a  block  diagram  of  a 
generic  LQG  controller.  In  the  context  of  the  adaptive  optics  control  of  this  research,  the 
z(tj-)  are  the  discrete-time  measurements  ftom  a  Hartmann  slope  sensor,  and  the  u (it)  are 
the  mirror  control  voltages. 


KALMAN  FILTER 


Z(i.) 


t  fc/^V—  J  Kit  ) 

l  _ ^  > _  .  .  1  j 

1  ►vJ 

1  H 
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— .  ^(<r  m  + 
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mor  j 
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i 

Figure  2.7.  Block  Diagram  of  LQG  Controller 

2.4.I  Kalman  Filter  The  Kalman  filter  for  the  simple  LQG  controller  can  be  de¬ 
signed  separately  from  the  regulator  due  to  the  certainty  equivalence  principle  [25:17], 
The  Kalman  filter  accepts  as  inputs  the  measurements  z(t,)  and  generates  an  estimate  of 
the  system  state,  xft,).  It  accomplishes  this  by  using  an  internal  model  of  the  system  to 
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propagate  this  state  estimate  and  its  covariance,  P  ,  and  then  updating  its  estimate  by 
appropriate  processing  of  measurements.  The  filter  covariance  is: 


P  =  £{[x-ic][x-xjT} 


(2.42) 


As  will  be  seen  shortly,  in  a  Kalman  filte:  this  covariance  is  not  dependent  on  the  actual 
values  of  the  measurements. 

Kalman  filter  operation  can  be  divided  into  two  sequential  processes  involving  the 
filter  state  estimate  and  covariance:  propagation  and  measurement  update.  Propagation 
is  the  change  in  the  state  estimate  and  filter  covariance  as  measurement-free  time  elapses. 
The  measurement  update  process  is  the  change  in  the  state  estimate  and  filter  covariance 
as  a  set  of  measurements  is  processed.  The  governing  equations  are  [25:18-19]: 


Hh  ;  = 

${U 

~  U-})ik{it-i)  t  ™dM.{U-i) 

(2.4.H) 

p  (tn  = 

Hu 

-  -  U- 1)  +  GdQdGj 

(2.44) 

K(*0  = 

p(*r 

)ht  [hp(<~)ht  +  r]  ~l 

(2.45) 

wi)  = 

Ut? 

)  +  K(t,)  [z,  -  Hx(f")] 

(2.46) 

p  (tf)  = 

p(t- 

■)-K(<i)HP(<r) 

(2.47) 

where 


m ) 
p  (tn 

K  (U) 

A 

Utf) 

m) 


estimate  of  x  at  start  of  t,-_ i  — »  t,  propagation  cycle 
filter  covariance  at  start  of  r,_]  — *  t,  propagation  cycle 
estimate  of  x  at  end  of  propagation  cycle  at  time  t, 
filter  covariance  at  end  of  propagation  cycle  at  time  t , 
Kalman  filter  gain  at  update  time 
measurement  realization  at  update  time 
estimate  of  x  just  after  update 
filter  covariance  just  after  update 


The  Kalman  filter  gain  equation  implies  an  on-line  inversion  of  an  m  x  m  matrix. 
If  the  number  of  measurements  is  greater  than  the  number  of  states  (m  >  n)  then  an 


2-20 


alternative  measurement  update  may  be  preferred  [26:54]: 


P(r  )  =  [p_1(/“)  +  HTR_1Ilj  ’  (2.48) 

K(tt)  =  P(^)HTR-1  (2.49) 

This  form  implies  on-line  inversions  of  matrices  having  dimensions  of  n  x  it. 

2.4.2  Linear  Quadratic  Regulation  As  derived  by  Maybeck  [25]  using  a  dynamic 
programming  approach,  the  optimal  control  input  can  be  written  as: 

H*  (*(«+),<,-)  =  -G*e(ii)Uit)  (2.50) 

Assuming  time-invariant  cost  matrices,  the  optimal  controller  gain  G*(t,)  is  : 

G;(/,)=  (u  +  BjKc(t,+i)Bj] _1  BjKc(t,+i)$(fi+i  -  U)  (2.51) 

where  the  controller  Riccati  matrix  Kc(t,  )  is  propagated  via  the  backward  Riccati  difference 
equation  [25:15]: 


K r(i.)  =  X  +  <f.r(<lTl  -  u) Kc(*i+1)*(li+1  -  U) 

-*T(*,+i  -  ti)Ke{U^)Bd  [U  +  BjKc(t1+] )Bd] _1 
xBjKc(i1+1)4*(t,+1  —  i,-)  (2.52) 


from  the  terminal  condition: 

Kc(<a+i)  =  X/  (2.53) 

This  backward  Riccati  equation  generally  exhibits  a  terminal  transient  as  tjv+i  is  ap¬ 
proached.  If  the  system  is  never  expected  to  reach  the  terminal  condition,  i.e.,  if  tyv  +  ]  =  00, 
then  the  steady-state  solution  to  the  backward  Riccati  equation  can  be  used  for  Kc  for  all 
bounded  time.  Substituting  this  Kc  into  Equation  (2.51)  yields  the  steady-state  controller 
gain  G‘. 
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2.5  Summary 


Light  traveling  through  atmospheric  turbulence  is  subject  to  spatial  and  temporal 
variations  in  index  of  refraction  along  the  propagation  path.  These  variations  result  in 
phase  distortions  in  formed  images.  The  statistics  of  such  distortions  are  dependent  on 
temporal  and  spatial  characteristics  of  the  atmosphere.  The  structure  function  expresses 
the  random  nature  of  some  of  these  statistics,  and  is  the  foundation  upon  which  much  of 
the  literature  builds. 

A  quantitative  description  of  the  phase  distortion  in  an  image  can  be  expressed  as 
an  infinite  sum  of  weighted  Zernike  basis  functions.  That  the  phase  distortion  is  not 
necessarily  time  invariant  is  expressed  in  the  time-dependency  of  the  Zernike  coefficients. 

Linear  Quadratic  Gaussian  (LQG)  control  is  one  approach  to  controlling  a  dynamic 
system.  Such  a  system  is  expressed  in  state-space  form  as  a  set  of  linear  dynamics  equations 
and  a  set  of  linear  measurement  equations.  LQG  control  utilizes  estimates  of  the  system 
state  to  derive  a  set  of  control  inputs,  the  goal  being  minimization  of  a  cost  function.  For 
the  adaptive  optics  system,  the  pi.ase  distortion  in  an  image  reflected  from  a  deformable 
mirror  is  to  be  minimized. 
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III.  Stochastic  Models 


3.1  Adaptive  Optics  System 

The  adaptive  optics  system  of  interest  is  a  ground-based  telescope  which  is  used  to 
observe  artificial  earth  satellites.  It  is  assumed  that  a  monochromatic,  coherent  source  of 
514-nm  light  is  within  the  same  isoplanatic  patch.  It  is  also  assumed  that  the  light  intensity 
is  of  sufficient  strength  to  be  seen  against  background  light.  A  simplified  schematic  of  the 
adaptive  optics  system  is  shown  in  Figure  3.1.  The  key  components  to  be  modeled  and/or 
designed  are: 

1.  Atmospheric  distortion 

2.  Deformable  mirror 

3.  Wavefront  sensor 

4.  Kalman  filter 

5.  Controller 

The  development  which  follows  treats  these  topics  in  the  order  indicated.  The  overall 
system  states  will  consist  of  the  atmospheric  distortion  Zernike  coefficients  augmented 
with  the  Zernike  coefficients  corresponding  to  the  mirror  shape. 

3. 2  Atmospheric  Effects 

Based  on  Taylor’s  frozen  turbulence  concept  introduced  in  Section  2.1.3,  the  temporal 
statistics  of  image  distortion  can  be  modeled  as  the  spatial  statistics  “blowing  by.”  A 
common  analytical  way  of  implementing  these  temporal  statistics  is  to  let  the  phase  front 
distortion  in  the  receiving  aperture  be  modeled  using  Zernike  basis  functions,  and  let  the 
Zernike  coefficients  (i.e.  the  elements  of  where  the  subscript  a  denotes  atmosphere) 
be  outputs  of  shaping  filters.  [16].  According  to  Noll  [29:210],  these  Zernike  coefficients 
are  well-modeled  as  zero-mean  Gaussian  random  processes.  Gaussianness  results  from 
the  summation  of  the  distortions  from  each  atmospheric  layer  of  turbulence  traversed. 
The  turbulence  of  each  layer  contributes  a  random  increment  to  the  final  set  cf  Zernike 


ATMOSPHERIC  BEAMSPLITTER  imaging 


Figure  3.1.  Simple  Schematic  of  Adaptive  Optics  LQG  Control 

coefficients  of  phase  distortion.  The  central  limit  theorem  states  that,  the  sum  of  many 
such  independent  random  contributions  is  Gaussian. 

The  general  equations  describing  the  atmospheric  distortion  shaping  filters  are  of  the 

form: 

*a(0  =  FaX^t)  +  Gawa(f)  (3.1) 

The  Fq,  Ga,  and  Qa  are  modeled  in  this  research  as  time  invariant.  Absence  of  a  B0u(t) 
term  is  due  to  the  atmospheric  distortion  being  uncontrollable.  The  dimension  of 
will  in  general  be  greater  than  the  number  of  Zernike  modes  modeled,  since  each  Zernike 
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coefficient  may  require  multiple  shaping  filter  states. 

For  this  research  atmospheric  distortion  is  modeled  as  consisting  of  Zernike  inodes 
1-14.  The  processes  a^(t)  corresponding  to  these  14  coefficients  can  be  extracted  from 
the  total  distortion  state  vector  x^(f)  via  multiplication  by  an  extraction  matrix  A  : 

a«(0  =  Axjt)  (3.2) 

The  vector  a ^(t)  denotes  the  Zernike  coefficients  of  the  image  phase  distortion  (excluding 
piston)  entering  the  adaptive  optics  system. 

One  means  of  designing  atmospheric  distortion  shaping  filters  (i.e.  FQ,  Ga  and 
Qa)  would  be  to  plot  actual  power  spectral  density  (PSD)  fata  for  each  Zernike  coeffi¬ 
cient,  based  on  collected  data.  For  added  authenticity,  the  data  could  be  collected  using 
ground-based  sensors  at  the  actual  telescope  site.  Having  the  PSDs  plotted,  straight-line 
approximations  could  be  drawn  to  determine  corner  frequencies  and  strengths  of  shaping 
filter  driving  noise.  Actual  PSD  data  is  unavailable  for  this  research,  but  does  exist  [23], 
Glasson  and  Guha  [10:13]  use  this  approach  to  model  the  phase  distortion  due  to  atmo¬ 
spheric  turbulence  for  the  first  five  Zernike  modes.  They  accomplish  this  modeling  by 
fitting  a  series  of  straight  lines  to  simulated  power  spectral  density  data. 

Another  approach  for  designing  the  atmospheric  distortion  shaping  filters  is  to  sim¬ 
ulate  the  Zernike  coefficients’  autocorrelation  kernels  from  analytically-derived  equations, 
followed  by  curve-fitting  to  standard  shaping-filter  equations.  Derivation  of  the  required 
analytical  relations  is  beyond  the  scope  of  this  research,  as  Section  2.1.4  points  out.  Ap¬ 
pendix  D  presents  the  details  of  generating  the  simulated  autocorrelation  kernel  data  and 
the  use  of  curve-fitting  the  data  to  shaping  filter  functions.  Also  presented  in  Appendix  D 
are  the  matrices  F0,  and  Qq  of  Equation  (3.1)  as  well  as  the  extraction  matrix  A  of  Equa¬ 
tion  (3.2).  The  software  developed  to  generate  simulated  autocorrelation  data  is  archived 
at  A  FIT  [39], 


3.3  Mirror 


This  research  assumes  the  mirror  is  a  monolithic  deformable  mirror  with  129  (97 
active)  evenly-spaced  actuators,  manufactured  by  ltek.  A  sample  of  this  mirror  is  at  the 
optics  laboratory  at  the  AFWL.  The  mirror  consists  of  four  major  components:  faceshcet. 
base,  electronic  circuitry,  and  actuators.  The  facesheet  is  a  monolithic  piece  of  ultra-low 
exnansion  glass  (ULE).  The  129  actuatoi  pusher  pads  are  precisely  machined  into  the  back 
side  of  the  facesheet.  The  base  is  made  of  similar  ULE  to  minimize  any  relative  thermal 
expansion  ellccts,  i.e.,  lateral  forces  on  the  actuators  [21:1]. 

Each  of  the  129  piezoelectric  actuators  is  constructed  from  layers  of  lead  magnesium 
niobate  (PMN).  Each  actuator  is  epoxied  to  both  the  facesheet  and  the  base.  The  actuators 
are  electrorestrictive,  meaning  cither  a  positive  or  negative  voltage  causes  the  actuator  to 
contract,  i.e.,  the  piezoelectric  stack  to  “shorten.”  To  make  two-way  excursions  of  the 
acLuatoi  possible,  the  stack  is  hiased  to  -150  volts.  Thus,  application  of  a  positive  voltage 
reduces  the  total  voltage  magnitude,  causing  the  stack  to  expand.  Application  of  a  negative 
voltage  increases  the  total  voltage  magnitude,  causing  the  stack  to  contract.  The  maximum 
magnitude  of  applied  voltage  allowed  is  300  volts  [21:1]. 

Application  of  these  large  voltages  to  the  actuators  is  controlled  by  a  control  voltage 
whose  range  is  ±10  volts.  A  -10  volt  control  voltage  corresponds  to  -150  volts  applied; 
making  the  total  voltage  magnitude  300  volts.  A  +10  volt  control  voltage  corresponds  to 
+  150  volts  applied,  making  the  total  voltage  magnitude  0  volts.  Thus,  positive  control 
voltage  causes  the  stack  to  expand;  negative  control  voltage  causes  it  to  contract  [21:1]. 

Of  the  129  actuators  on  the  mirror,  only  the  central  97  are  actively  controllable. 
The  ltek  Operation  Manual  [21:4]  states  the  remaining  32  could  be  made  independently 
controllable  by  the  addition  of  respective  driver  electronics.  At  present,  these  extras  are 
tied  to  a  bias  voltage  to  provide  fixed  boundary  conditions.  Figure  3.2  shows  the  location 
and  numbering  scheme  for  the  97  active  actuators,  assuming  the  manufacturer-specified 
0.85-cm  spacing. 

It  should  be  noted  that  the  configuration  of  this  mirror  at  the  AFWL  optics  develop¬ 
ment  laboratory  has  only  the  central  09  actuators  independently  controlled;  the  remaining 
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Figure  3.2.  Actuator  Locations  on  Itek  S7-Actuator  Mirror 


28  are  slaved  to  the  nearest  neighbor  within  the  central  69  [23].  This  arrangement  is  used 
to  implement  the  so-called  “zonal”  approach  to  mirror  control  in  which  each  actuator  cor¬ 
responds  to  a  spatial  sample  of  the  image  in  the  measurement  device.  In  this  case,  the 
Hartmann  wavefront  sensor  (to  be  discussed  in  Section  3.5)  has  69  subapertures.  This 
research,  however,  uses  the  so-called  “modal”  approach,  which  implies  control  of  (Zernike) 
modes  of  the  distortion.  Greenwood  [11:549]  states,  without  justification,  that  the  de¬ 
grees  of  correction  possible  with  either  approach  are  similar.  He  further  states  the  current 
(197S)  preference  is  the  zonal  approach.  Southwell  [35:1006],  on  the  other  hand,  argues 
that  reconstruction  of  the  phasefront  from  Hartmann-type  sensor  measurements  appears 
to  be  superior  for  the  modal  approach.  This  research  will  not  compare  one  approach  to 
the  other;  the  zonai  approach  is  merely  mentioned  to  explain  the  AFWL  configuration. 
This  research  assumes  all  97  active  actuators  are  driven. 
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Tiie  mirrored  surface  of  the  facesheet  in  front  of  the  97  actuators  lies  within  a  4.25- 
cm  radius  opening  in  a  circular  bezel.  Tims  tlu*  effective  diameter  of  the  mirror  aperture 
is  8.5  cm.  Table  3.1,  condensed  from  the  It.ek  Operation  Manual  [21:6],  shows  add  'ional 
mirror  parameters. 


Table  3.1.  AFWL  Deformable  Mirror  Characteristics 


Faccsheot  material 

ULE 

Clear  aperture 

8.5  cm 

Number  of  actuators 

97  controlled  +  32  biased 

Actuator  geometry 

Square  array,  lr  across  diameter 

Actuator  spacing 

0.85  cm 

Hysteresis 

None  observed 

Stroke 

3.91  microns  (mean) 

Surface  figure 

A/10  p-p  @  0.6328  fi 

Coating 

Protected  aluminum 

Reflectivity 

>84% 

Actuator  bandwidth 

500  Rz 

Operating  temperature 

20°C-30°C 

Package  size 

12-inch  cube 

Linearity  of  a  single  actuator  implies  that  the  graph  of  local  mirror  displacement 
versus  control  voltage  is  a  straight  line.  AFWL  performed  such  measurements  on  a  few 
actuators;  one  of  their  plots  is  shown  in  Figure  3.3  [22].  Because  of  the  flattening  at  the 
ends  of  the  control  voltage  range,  the  behavior  is  not  strictly  linear.  However,  over  a 
reasonable  range  of  control  voltage,  the  function  is  approximately  linear.  The  design  of 
the  controller  assumes  a  linear  model.  It  is  expected  that  an  appropriate  weighting  matrix 
in  the  cost  function  will  keep  the  control  voltage  within  the  linear  range  most  of  the  time. 

Linearity  of  the  colicctionoi  actuators  implies  inter-actuator  superposition  holds.  For 
example,  applying  one  control  volt  to  actuator  #  39,  measuring  the  mirror  response,  then 
applying  one  control  volt  to  actuator  #  59.  and  adding  the  responses  should  yield  the  same 
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Figure  3.3.  Sketch  of  Actuator  Displacement  Linearity 

total  response  as  applying  the  voltages  simultaneously.  Such  multi-actuator  measurement 
data  are  not  available  for  this  mirror,  but  this  research  will  assume  superposition  holds. 
The  literature  contains  precedent  for  the  validity  of  the  superposition  assumption  [1:31]. 

3.3.1  Steady-State  Mirror  Behavior  Ii  is  expected  and  desired  that  application  of 
a  control  voltage  to  an  actuator  will  cause  the  mirror  to  deform.  This  deformation,  How¬ 
ever,  does  not  occur  instantaneously;  a  finite  time  is  required  for  the  mirror  to  attain 
its  “steady-state”  position.  This  section  of  the  report  addresses  the  relationship  between 
applied  control  voltage  and  the  steady-state  mirror  response.  The  eventual  goal  is  to  re¬ 
late  a  vector  of  97  control  voltages  to  a  vector  of  14  resulting  Zernike  coefficients.  These 
“mirror”  Zernike  coefficients  add  to  the  14  “atmospheric”  Zernike  coefficients  of  the  in¬ 
cident  atmosphere-aberrated  light,  hopefully  cancelling  them.  Thus,  reflection  from  the 
deformable  mirror  is  modeled  as  the  addition  of  atmospheric  Zernike  coefficients  to  the 
mirror  Zernike  coefficients. 


3-7 


The  influence  function  is  a  mathematical  representation  of  the  effect  of  a  single  ac¬ 
tuator  voltage  on  the  local  mirror  shape.  Usually,  the  influence  function  is  nonzero  only 
in  the  vicinity  of  the  actuator;  the  influence  function  of  an  actuator  has  a  limited  spa¬ 
tial  domain.  Several  factors  affect  the  influence  function  of  an  actuator.  These  factors 
include  piezoelectric  type,  facesheet  material  and  thickness,  proximity  and  geometry  of 
neighboring  actuators,  and  actuator  linearity.  Based  on  the  previous  assumptions  of  actu¬ 
ator  linearity  and  evenly-spaced  geometry,  it  is  assumed  that  all  actuators  have  the  same 
influence  function.  The  fact  that  the  word  “mean”  appears  in  Table  3.1  for  the  actuator 
stroke  parameter  implies  that  each  actuator  does  not  have  the  same  influence  function. 
Nevertheless,  this  research  assumes  uniformity  of  the  97  influence  functions. 

Limited  influence  data  is  available  for  the  Itek  mirror  at  the  AFWL.  Itek  provided 
data  corresponding  to  a  two-dimensional  slice  of  the  influence  function  for  the  central 
actuator,  #49.  This  slice  was  taken  along  the  mirror  X-axis,  spanning  actuators  71,  60, 
49,  38,  and  27.  Itek’s  plot  of  this  data  is  shown  in  Figure  3.4.  The  magnitude  of  the  voltage 
applied  to  the  actuator  was  200  volts  for  this  data.  Influence  data  was  also  available  along  a 
45-degree  slice  through  the  central  actuator.  This  data  showed  less  of  a  negative  excursion, 
but  was  of  the  same  general  shape  as  the  slice  along  the  X-axis.  Thus,  it  is  assumed  that 
the  slice  through  the  X-axis  is  approximately  valid  in  any  direction.  The  ordinate  data 
(X)  used  to  plot  Figure  3.4  were  converted  into  units  of  cm.  The  abscissa  data  were 
converted  into  2-way  decrease  in  path  length,  measured  in  wavelengths  of  514-nrn  light. 
The  abscissa  data  were  then  scaled  to  represent  one  volt  of  control  voltage,  invoking  the 
linearity  assumption.  The  transformed  data  were  curve-fit  to  a  tenth-order  polynomial  in  X 
truncated  at  X=±  1.6739  cm,  the  spatial  domain  of  influence.  The  truncated  function  was 
then  rotated  about  the  actuator  axis  (i.e.,  the  “Z”  axis)  resulting  in  the  three  dimensional 
influence  function  of  the  central  actuator.  This  influence  function,  when  expressed  in 
mirror  coordinates,  yields; 

f(X,Y\XA,YA)  =  {0.9673  (3.3) 

-2.726  [(X  -  AY)2  +  {Y  -  A*)2] 

+2.943  [(A  -  A'a)2  +  (>'  -  F/t)2]2 
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Position  (inches  --  resilngs  taken  every  0.020') 

Figure  3.4.  X-Axis  Slice  of  Actuator  #49  Influence  Function 

-1.573  [(A'  -  XV7  t  (F  -  Fa)2]3 
+0.4137  [(A-  X„)2  +  (K-y,02]4 
-0.04236  [(A  -  AA)2  +  (y-yA)2j5 } 
x  ([1.6739]'  -  [(A  -  A Af  +  (F  -  F*)2])}  (3. 

where 

A,  F  =  coordinates  of  point  on  mirror  (cm) 

Xa,  Ya  =  coordinates  of  actuator  (cm) 

u._i(-)  =  unit  step  function 

/  =  2-way  path-length  decrease  (wavelength/ volt) 
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Figure  3.5  shows  the  three-dimensional  piot  of  the  influence  function  for  the  central  actu¬ 
ator.  Note  the  X,Y  units  are  cm,  whereas  the  “Z”  units  are  wavelengths.  This  accounts 
for  the  extreme  protruding  appearance  of  the  plot. 


Figure  3.5.  Approximate  Influence  Function  for  Actuator  #49 

At  this  point,  a  mental  example  is  helpful.  Suppose  a  perfectly  pla.nar  wavefront  is 
incident  on  the  mirror.  Further,  suppose  that  one  volt  of  control  voltage  is  applied  to  an 
arbitrary  actuator,  ail  other  control  voltages  being  zero.  Let  the  X,Y  position  of  interest 
be  right  at  the  actuator,  i.e.,  X  —  X 'a  and  Y  =  YA.  Therefore,  the  value  of  /  at  that 
actuator  is  0.9673  wavelengths.  Actually,  the  mirror  only  moves  0.4836  wavelengths  at  the 
actuator,  but  this  causes  the  path-length  of  the  reflected  light  to  be  reduced  by  twice  as 
much  and  /  is  defined  to  be  such  reduction  in  path  length.  Suppose  the  light  reflected 
from  the  entire  mirror  then  passes  through  an  aperture  of  the  same  "ameter.  The  light 
in  this  aperture  corresponding  to  that  reflected  from  the  extended  actuator  reaches  the 
aperture  first  since  its  path  length  is  shorter.  Thus,  its  absolute  phase  is  a  larger  value 
than  that  of  the  liyht  from  the  rest  of  the  miiror.  Therefore,  a  plot  of  the  phase  deviation 
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in  the  aperture  looks  exactly  like  the  influence  function  plot  with  piston  subtracted  out. 

The  previous  discussion  can  be  extended  to  the  case  where  control  voltages  are  ap- 
plii  -1  to  all  97  actuators.  Assume  that  tl  e  light  incident  on  the  mirror  is  perfectly  planar. 
Thus  any  phase  distortion  in  the  reflected  image  will  be  solely  due  to  the  mirror  shape. 
All  97  actuators  are  then  excited  by  arbitrary  but  time-invariant  control  voltages.  When 
the  mirror’s  shape  reaches  steady  state,  a  snapshot  of  the  phase  distortion  in  the  reflected 
image  can  be  obtained  by  subtracting  the  spatial-average  phase  (i.e.,  the  piston).  An 
appropriate  equation  for  the  phase  distortion  in  the  reflected  image  is: 


4>(X,Y)  = 


97 


£  f(X,Y.XMYA)u(A) 


L/l  =  l 


-£ 


97 


£  !{X,y,xa,ya)u{a) 


>4=1 


(3.5) 


Equation  (3.5)  represents  the  mirror’s  contribution  to  the  distortion  of  the  reflected  image. 
Since  the  left-hand  side  is  phase  deviation  in  an  aperture,  it  can  be  expressed  as  a  linear 
combination  of  Zernike  functions  Zx  through  Z <*>.  For  this  research,  the  atmospheric 
distortion  is  modeled  as  the  first  14  of  these.  It  then  follows  the  mirror  will  only  be  trying 
to  correct  only  the  first  14  modes.  Thus  the  significant  phase  deviation  caused  by  some 
operational  set  of  commands  to  the  mirror  car.  be  modeled  as  linear  combinations  of  the 
first  14  Zernike  functions  (sis  usual,  excluding  piston): 


<*>(*,  y)«£  a,  Zi{X,Y) 
!=.! 


(3.G) 


The  a,  Zernike  coefficients  can  be  determined  using  Equation  (2.25),  expressed  below  in 
rectangular  form: 

„  _  /  dY  rdXW(X,Y)<p{X,Y)Z,(X^Y) 

JdYfdXW(X,Y )  1  j 

Substituting  Equation  (3.5)  into  Equation  (3.7),  and  realizing  that  the  product  of  piston 

and  any  non-piston  Zernike  function  integrates  to  zero  yields: 


I  dY  f  dXW(X,Y )  \ZA=,  f(X,Y,XA,YAHA)\  Z{(X,Y ) 
J  dY~f  d X  W(X,Y) 


(3.8) 
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One  can  then  then  pull  the  summation  outside  the  integrals,  and  commute  the  scalar 
control  voltage: 


«<  =  £{ 

^=1  k 


97  dY  f  dXW(X,Y)f(X,Y,XA,YA)Zi(X,Y)l 


J  dY  I  dX  W(X,Y) 
Equation  (3.9)  can  be  expressed  in  vector  form: 

T 

a.  =  m,  u 


(3-9) 


(3.10) 


where  the  A-th  comp<  nent  of  the  m  vector  is  the  projection  of  the  A-th  actuator’s  influence 
function  along  the  i-th  Zernike  function.  Likewise  the  A-th  component  of  the  u  vector  is 
the  control  voltage  on  the  A-th  actuator.  Equation  (3.10)  can  be  written  for  each  Zernike 
function: 


T 

fll 

=  Eli  M 

0-2 

SI 

SI 

II  •• 

«14 

T 

=  mu  n. 

(3.11) 


Finally,  Equations  (3.12)  can  be  combined  into  matrix  form  as: 


a  =  Mu 


(3.12) 


where  a  is  the  vector  of  Zernike  coefficients  describing  the  mirror’s  steady-state  contribution 
to  the  distortion  of  the  reflected  image.  The  matrix  M  is  the  steady-state  influence  matrix, 
since  it  relates  the  steady-state  influence  of  the  mirror  to  the  applied  control  volt  ages.  Each 
element  of  the  matrix  is  the  projection  of  an  actuator’s  influence  function  along  a  Zernike 
function  direction.  The  Fortran  program  which  calculates  the  M  matrix  for  this  research 
is  archived  at  AFIT  [39].  The  resulting  M  matrix  is  shown  in  Appendix  E. 

In  order  to  check  the  reasonableness  of  the  steady-state  influence  matrix,  it  was 
decided  to  perform  a  test.  The  test  consisted  of  analytically  determining  the  required 


3-12 


control  voltage  vector  to  cause  the  reflected  image  to  consist  of  a  selected  Zernike  mode, 
assuming  perfectly  planar  incident  light.  Then  the  resulting  voltage  vector  was  applied 
to  the  mirror  in  a  mathematical  simulation  to  see  if  the  resulting  phase  deviation  plot 
“looked  like”  a  plot  of  the  selected  Zernike  function.  The  14-th  Zernike  function  was 
selected  because  of  its  high  spatial  frequency  content  relative  to  the  other  modes.  If  the 
mirror  is  able  to  reproduce  the  14-th  mode  well,  it  should  be  able  to  reproduce  lower-  order 
modes  as  well.  An  arbitrary  value  of  fli4  was  chosen,  with  the  remaining  coefficients  set  to 
zero.  Equation  (3.12)  was  then  solved  for  the  required  voltage  vector[l4:35]: 

u  =  MT(MMT)'1a  (3.13) 

This  is  the  unweighted,  rniniinum-norm  solution.  This  series  of  matrix  and  vector  op¬ 
erations  resulted  in  a  set  of  97  control  voltages.  A  Fortran  program  was  then  used  to 
simulate  the  phase  deviation  caused  by  the  set  of  voltages  and  generate  plotable  results. 
Comparison  of  Figure  3.6  with  the  plot  of  the  14-th  Zernike  function  in  Appendix  A  lends 
credibility  to  the  calculated  M  matrix. 


Figure  3.6.  Simulated  Mirror  Reconstruction  of  14-th  Zernike  Mode  Phase  Distortion 


3-3.2  Transient  Behavior  Up  until  now,  only  the  steady-state  behavior  of  the  mir¬ 
ror  has  been  considered;  the  mirror  has  been  given  time  to  respond  completely  to  the  input 
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control  voltages.  Now  the  transient  behavior  of  the  mirror  wili  be  considered.  The  mirror's 
manufacturer  description  states  that  an  actuator  acts  as  a  electrical  capacitive  load  [21]. 
It  is  thus  reasonable  to  expect  that  application  of  a  step  control  voltage  to  an  actuator 
does  not  result  in  a  step  response,  but  rather,  an  exponential  approach  to  an  asymptote. 

The  Operation  Manual  [21:6]  for  the  Itek  deformable  mirror  indicates  the  mean  band¬ 
width  for  the  actuators  is  500  Hz.  AFWL  bandwidth  data  [22]  is  generally  in  agreement. 
Figure  3.7  is  the  AFWL  plot  of  the  frequency  response  for  one  of  the  97  active  channels. 
For  this  particular  actuator  channel,  the  -3-dB  frequency  is  497  Hz.  It  is  assumed  this 


X-497.  2  Ha 
V„--3.  (307  1  d9 
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V  — —  3.  0007  dB 


Figure  3.7.  Frequency  Response  of  An  AFWL  Mirror  Actuator  Channel 


plot  includes  the  effects  of  the  driver  circuitry  dynamics  and  the  mirror  dynamics.  The 
approximate  slope  of  the  rolloff  is  -10  dB  per  decade  of  frequency,  which  confirms  the 
presence  of  a  first  order  pole.  It  is  assumed  the  decibel  values  plotted  are  normalized  with 
respect  to  the  low  frequency  response: 


#  of  dB  --  10  log 


Mjf 


(3.14) 
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where  R  is  the  the  ratio  of  sinusoidal  excitation  magnitude  to  sinusoidal  response  magni¬ 
tude.  Assuming  the  dynamics  of  the  mirror  are  first  order,  time-invariant,  and  determin¬ 
istic,  the  dynamics  equation  for  the  displacement  at  the  actuator  site  tor  control-voltage 
excitation  of  the  actuator  can  be  modeled  by  the  scalar  dynamics  equation: 


g(t)  =  c  g(t)  +  d  u(l) 


(3.15) 


where  g  is  the  decrease  in  2-way  path-length  of  light  reflected  at  the  actuator  (i.e.,  twice 
the  actual  physical  displacement  of  the  mirror  at  the  actuator  site).  The  values  of  c  and  d 
in  Equation  (3.15)  will  now  be  analytically  determined,  using  the  bandwidth  value  and  the 
peak  displacement  of  the  assumed  influence  function,  Equation  (3.4).  Taking  the  Laplace 
transform  of  Equation  (3.15)  for  zero  initial  displacement  and  rearranging  yields: 


G'(5)  = 


(3.16) 


The  s  is  the  Laplace  complex  variable.  Now  if  the  control  voltage  input  is  switched  from 
zero  to  a  sinusoidal  voltage  at  time  t=0,  the  input  Laplace  transform  can  be  found  from  a 
one-sided  Laplace  transform  table  [8:772] : 

u(t)  —  sinwf,  t>0  (3-17) 

U(S)  =  (3’18) 

s'  +  u r 


Substituting  this  into  Equation  (3.16)  and  again  referring  to  a  Laplace  transform  table 
[8:772]  the  time-domain  response  can  be  obtained: 


+ 


cl  i  u'  wV^-fw2 


:  sin 


^ ut  — 


tan 


-l 


w 

-c 


,  t  >  0 


(3.19) 


Since  c  <  0,  the  exponential  term  in  the  above  equation  vanishes  as  t  -*  oo.  Thus,  after 
the  transient  response  has  died  out,  the  sinusoidal  “steady-state”  response  is: 


<?«(*)  =  —  =f==  sin  ( wt  -  tan-1  [—  ^  ,  t  >  0 

VC'  t  UJ~  \  l  ~  Cj  / 


(3.20) 
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At  zero  frequency,  the  amplitude  of  g  is  merely  At  the  bandwidth  frequency  (assumed  to 
be  500  Hz,  or  10007T  rad/sec),  the  amplitude  of  g  should  be  Therefore  the  equation 


J_  d  _ _ d _ 

y/2  c~  y/c*  +  (lOOOff)2 


(3.21) 


can  b  solved  for  c,  yielding  : 


IOOOtt 


c  —  — 


sec 


-1 


y/2 
-2222  sec-1 
1 

0.00045  sec 


(3.22) 


Thus  it  can  be  said  the  time  constant  of  a  typical  actuator  is  approximately  0.45  millisec¬ 
onds. 

From  Equation  (3.4)  the  steady-state  value  of  g  for  a  unit  step  control  voltage  u(t) 
is  0.S673  wavelengths.  The  left-hand  side  of  Equation  (3.15)  is  zero  at  steady-state  for  a 
step  input.  With  zero  on  the  left-hand  side  and  uss(t)  set  equal  to  one,  Equation  (3.15) 
can  be  solved  for  d.  The  resulting  value  for  d  is: 


««(*) 
-(- 


'(_17r)  (°-9673) 


Thus  the  dynamics  of  an  actuator  can  be  approximated  by: 

g(t)  = - - - g(t)  +  2150  u(t) 

v  ’  0.00045  JW  w 


(3.23) 


The  validity  of  Equation  (3.23)  for  any  actuator  obviously  depends  on  the  assumptions 
of  all  actuators  having  the  same  influence  function  and  bandwidth.  Also  the  dynamics 
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are  assumed  to  be  deterministic.  With  additional  statistical  data  on  the  bandwidths  and 
influence  functions  of  all  actuators,  97  versions  of  Equation  (3.23)  could  be  written.  Also, 
the  noise  in  th<‘  control  voltage  could  b>  modeled,  causing  the  addition  of  a  stochastic  term 
to  the  equations.  This  research  assumes  the  deterministic  Equation  (3.23)  is  valid  for  all 
actuators. 

The  previous  discussion  of  dynamics  at  the  actuator  site  can  be  extended  to  the  phase 
distortion  over  the  entire  mirror.  It  is  desired  to  formulate  the  dynamics  of  the  mirror  in 
the  form  of: 

=  f’mxm(t)  -r  Bmu(t)  (3.24) 

where  the  elements  of  vector  xm(t)  are  the  time-varying  Zernike  coefficients  which  de¬ 
scribe  the  mirror’s  contribution  to  the  phase  distortion  in  the  reflected  image.  The  lack 
of  stochastic  terms  is  an  engineering  assumption.  The  previous  assumption  of  identical 
actuator  bandwidths  implies  the  entire  mirror  has  the  same  bandwidth,  and  therefore  the 
same  time  constant.  If  a  set  of  97  control  voltages  is  simultaneously  applied  to  the  mirror, 
therefore,  the  dynamic  behavior  of  the  entire  mirror  surface  would  exhibit  behavior  indica¬ 
tive  of  a  0.00045-second  time  constant.  Thus,  if  <pm(A',Y,  <)  is  the  mirror’s  contribution 
to  the  phase  distortion  in  the  reflected  image  then  the  following  scalar  equation  can  be 
written: 

4>m(X,Y,t)  =  -  i  ^n(A\r,0+bI(A',r)u(0  (3.25) 

The  vector  of  functions  bm(A\  Y)  maps  the  voltage  of  each  actuator  to  the  rate- of- change 
of  phase  at  point  (X,Y).  The  phase  distortion  introduced  by  the  mirror,  $m(A',Y,f).  can 
be  approximated  by  a  linear  combination  of  a  finite  number  of  Zernike  basis  functions: 

4>m(X.  Y,  t)  ~  a^(t)Zi(X,  1 r)  -f  02(f)^2(A ,  Y)  + - f  oi4(0^m(A  ,  Y)  (3.26) 


Taking  the  time-derivative  of  Equation  (3.26)  yields: 


<£m(A ,  V,  t)  ~  di(t)Zi(A ,  Y)  t  d.2{t)Z-2(X ,  V  )  -f  -  •  •  +  dn(t)Zi4(A  ,  Y) 


(3.27) 


An  immediate  temptation  is  to  write  these  two  equations  in  vector  form,  followed  by 
substitution  into  Equation  (3.25),  which  would  yield: 

a(t)TZ(X,y)  =  -ia(0TZ(A',Y)  +  bl(X,Y)n(t)  (3.28) 


To  get  Equation  (3.28)  into  the  form  of  Equation  (3.21)  with  the  Zernikc  coefficients 
being  the  states,  one  is  further  tempted  to  take  the  transpose  of  both  side-  of  Equa¬ 
tion  (3.28),  then  premultiply  both  sides  by  the  vector  Z(A',y),  then  premultiply  both 
sides  by  the  matrix  inverse  [z(A, y)Z^(A',  V )j  .  The  problem  with  this  approach  is 
that  [z(A',y)ZT(A',y)l  is  a  rank-one  matrix  whose  inverse  does  not  exist. 

To  get  Equation  (3.25)  into  the  form  of  Equation  (3.24),  use  is  made  of  the  orthog¬ 
onality  property  of  the  Zernike  functions.  Substituting  Equations  (3.26)  and  (3.27)  into 
Equation  (3.25)  yields: 


X>(t)z,(A,y)  =  -i£>(ow,y)+ 1;  w*.i>j(o  (3-29) 

i=i  i=i  j= i 


Multiplication  of  Equation  (3.29)  by  Z\(X,Y),  then  spatially  integrating  both  sides  of  the 
equation  over  the  area  of  the  mirror  aperture — using  Equation  (2.22)  where  appropriate — 
yields: 

di(<)  =  --gj(0  +  (3.30) 

The  vector  £  is  defined  such  that  its  j-th  element  is  the  projection  of  the  function  bj( A',  Y) 
into  Zi(A'.y).  Multiplying  Equation  (3.29)  by  each  successive  Zernike  function  and  fol¬ 
lowing  the  integration  procedure  yields  the  set  of  equations: 


«i(0  =  -^ci(0  +  ^2/f7a(0 
O2(0  =  ~^a2{t)  +  ~/rJu(t) 

1  l  T‘ 

«14(  0  =  --om(0+ 


(3.31) 


(3.32) 


3- IS 


Putting  this  set  of  differential  equations  into  state-space  form,  yields: 


=  ding 


2W(0  +  Bmu(0 


(3.33) 


where  diag  [-]  is  a  diagonal  matrix  of  the  argument  and  the  rows  of  matrix  Bm  are  the  //^ 
row  vectors,  divided  by  7r R2. 

To  determine  the  elements  of  the  matrix  Bm,  consider  the  entire  mirror  at  steady 
state  such  that  xm(£)  =  0.  Substituting  this  into  Equation  (3.33)  and  manipulation  yields: 


=  diag[r]Bm«i„(£) 


(3.34) 


Comparing  this  equation  with  Equation  (3.12),  it  is  straightforward  to  show  that: 


B. 


-  diag 


II 

LrJ 


M 


(3.35) 


where  M  is  the  steady-state  influence  matrix  and  r  is  0.00045  seconds. 

In  summary,  the  dynamics  of  the  Itek  deformable  mirror  are  written  in  terms  of 
the  mirror’s  contribution  to  the  Zernike  coefficients  of  the  reflected  image.  Invoking  the 
assumptions  of  identical  actuator  influence  functions  and  identical,  first-order  actuator 
dynamics,  the  equation  for  these  contributions  to  the  Zernike  coefficients  is: 


xm(£)  =  diag 


11 

r 


xm(£)  +  diag 


Mu(() 


(3.36) 


It  should  be  noted  this  equation  is  a  truncation  of  the  inflnite-dirnensional  functional-space 
description  to  14  dimensions. 


3. J,  Augmented  System  Dynamics 

The  overall  continuous-time  description  of  the  adaptive  optics  system  dynamics  is: 


*a(0 

= 

Fa 

0 

£ 

4- 

0 

«(0  + 

3£a(  0 

_  £m(0  _ 

0 

Fm 

0 
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The  upper  partition  of  Equation  (3.37)  represents  atmospheric  dynamics;  the  lower  parti¬ 
tion  represents  mirror  dynamics. 


3.5  Wavefront  Sensor 

The  image  reflected  from  the  deformable  mirror  contains  phase  distortion  due  to 
atmospheric  turbulence  and  also  counterdistortion  duo  to  the  mirror  deformation.  Dis.or- 
tion  here  is  deviation  in  phase  from  the  aperture-average  phase,  thus  piston  is  excluded. 
Letting  the  phase  distortion  in  the  reflected  image  be  modeled  in  the  Zernike  functional 
space,  the  corresponding  discrete-time  Zernike  coefficients  for  the  reflected  image  are: 


«(U) 


A  |  I 


,(*<) 


(3.38) 


where  A  is  the  extraction  matrix  defined  in  Equation  (3.2)  and  I  is  a  14-by-14  identity 
matrix.  Note  that  this  equation  models  the  reflection  process  as  a  summation  of  two  sets 
of  Zernike  coefficients.  Since  the  wavefront  sensor  measures  the  distortion  in  the  reflected 
image,  the  discrete  measurement  equation  for  the  wavefront  sensor  can  be  written  as: 


H'a(t,)  +  v(N) 


H' 


A  |  I 


Hx(t.)  +  v(<.) 


*«(L)  I 


+  x(N) 


(3.39) 


This  research  takes  the  approach  of  first  calculating  the  H'  matrix  analytically,  then  post- 
multiplying  by  the  [  A  |  I  ]  matrix  to  obtain  H.  The  measurement  noise  covariance  R 
corresponding  to  v(t,)  is  left  as  a  design  parameter. 


3.5.1  Wavefront  Sensor  Description  The  reflected  image  is  diverted  to  the  wave- 
front  sensor  via  a  beamsplitter.  Light  reflected  from  the  beamsplitter  enters  the  Hartmann 
wavefront  sensor.  The  Hartmann  sensor  consists  of  an  array  of  square  convex  lenses,  each 
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X  (cm) 

Figure  3.8.  Subaperture  Locations  on  Hartmann  Wa'efront  Sensor 

considered  a  sub?  rture  of  the  sensor.  The  sense;  has  69  fully-illuminated  square  sub¬ 
apertures.  Figure  3.8  shows  the  arrangements  of  the  subapertures  in  the  Hartman  sensor 
aperture.  Each  sulaperture  is  0.0G  cm  on  a  side  [27]  and  focuses  its  share  of  the  incident 
light  onto  a  reticon  detector.  The  location  of  the  focused  spot  of  light  on  the  detector 
is  an  indication  of  the  average  x-  and  y-tilt  in  the  subaperture.  Since  subaperture  tilts 
are  measured,  the  Hartmann  sensor  is  essentially  a  slope  sensor.  Figure  3.9  shows  a  two- 
dimensional  analog  of  the  operation  of  the  Hartmann  sensor  [27]. 

3.5.2  Derivation  of  H'  An  important  concept  is  that  phase  distortion  in  the  image 
manifests  itself  as  a  set  of  subaperture  tilt  measurementss  in  the  Hartmann  sensor.  Thus 
any  Zernike  mode  (except  piston)  distortion  of  the  whole  image  is  measured  as  a  set  of 
tilt  measurements.  It  is  assumed  the  Hartmann  sensor  can  output  these  138  tilt  measure¬ 
ments  (69  x-tilts  and  69  y  -  tilts)  directly.  The  sensor  normally  processes  tie  set  of  138  tilt 
measurements  and  reconstructs  the  phase,  but  this  reconstruction  is  based  only  on  current 
tilts;  previous  time  realizations  are  not  considered  This  reconstruction  process  is  accom- 
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Figure  3.9.  Concept  of  Operation  for  Two-Dimensional  Hartmann  Sensor 

p'ished  in  software,  requiring  computational  time.  This  research  therefore  assumes  the  138 
noise-corrupted  tilt  measurements  are  the  outputs  of  the  sensor,  i.e.,  comprise  the  z(t,) 
vector.  The  Kalman  filter  estimates  the  system  states  from  the  138  till  measurements. 

Consider  the  phase  distortion  of  the  image  entering  the  sensor  to  be  a  function  of 
position  in  the  aperture  and  discretized  time: 

d>(A',  Y,  tj)  (3.40) 


The  variables  A'  and  Y  are  rectangular  coordinates  with  respect  to  the  center  of  the  total 
sensor  aperture.  Just  prior  to  light  entering  one  of  the  square  subaperture  lenslets,  the 
x-tilt  at  position  (A'i,Yj)  is: 


d(f>[X,Y,  tj) 
dY 


lA'j.y, 


(3.41) 


Letting  (A'j.V,)  define  the  coordinates  of  the  center  of  the  s-th  subaperture  with  respect 
to  the  entire  sensor  aperture,  the  average  x-tilt  going  into  the  s-th  subapertuie  lenslet  is 
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then: 


1  />'*+*?  /•*«  +  # 

-  /  ^  dY  *  dX 

A  J Yr-&  Jxt-& 


d<p{X,Y,  tj) 
dY 


where  A  is  the  area  of  the  square  subaperture. 


(3.42) 


According  to  Petersen  and  Cho  [32]  the  output  of  a  Hartmann-type  slope  sensor  is 
a  similar  integral,  but  with  the  slopes  spatially  weighted.  For  example,  the  output  of  the 
x-tilt  channel  in  a  subaperture  (ignoring  noise)  is: 


XI 


1 

^•laL 


.XI 

3 


dY 


-/  ««’.(* 


-  x„y-y,) 


D0(X,Y,t,) 

8Y 


(3.43) 


The  constant  L  is  dependent  on  wavelength,  lenslet  focal  length,  and  output  scaling.  Pe¬ 
terson  and  Cho  derive  the  spatial  weighting  function  \Vx(X',Y")  for  a  square  subaperture; 
the  normalized  version  is: 


Wx(X\Y')  = 


2 ln(2) 


"s ,  f  2y'hi 

to'1  +  73J  J 

(3-44) 


Similar  equations  can  be  written  for  the  subaperture  y-tilt  channel.  A  plot  of  Equa¬ 
tion  (3.44)  for  the  subaperture  dimensions  of  this  research  is  shown  in  Figure  3.10.  Obvi¬ 
ously  the  slope  of  light  traversing  the  center  of  the  subaperture  is  weighted  more  heavily 
than  that  traversing  near  the  edges. 

If  the  phase  distortion  of  an  image  can  be  represented  as  a  linear  combination  of 
basis  functions,  then  the  partial  derivative,  expressed  in  Equation  (3.41)  can  also: 


3<p(A',Y,  tj)  dZi(X,Y') 

'  -  =  Gl^J)  dY -  + 


a14(tj) 


az„(A\y) 

dY 


(3.45 


Substitution  of  Equation  (3.45)  into  Equation  (3.43)  and  simple  manipulation  yields: 


L 


14 


-(C)  =iE 


1=1 


fY-+^,y  /A*+^,v,f./v  v  v  v^z*x :, y ) 

°i y~.i i  I  —  d}  I  dX  If  i{-\  XS.Y  lj)'  „v 

JY.-'Sf-  J  x,-X±  dY 


(3.46) 
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Figure  3.i0.  Slone  Measurement  Weighting  Function 


Equation  (3.46),  representing  the  noise-free  output  of  the  x*tilt  channel  for  the  s-th  sub- 
aperture,  can  be  written  in  vector  notation: 


**.(*;)  =  ~  nj,  &(tj) 


(3.47) 


The  i-th  element  of  the  nI(  vector  can  be  thought  of  as  the  projection  of  the  i-th  Zernike 
function  into  x-tilt  measurement  space  at  the  s-th  subaperture.  Knowing  the  location  of 
the  center  of  an  actuator  relative  to  the  entire  aperture,  and  knowing  the  partial  derivatives 
of  the  Zernike  functions,  the  elements  of  n,.  can  be  determined  as: 

fY>  +  *r  fx>+ *r  dZiXY ) 

=  Jy^  ^  ^  W^X  -X^~  ^ 

Equation  (3.47)  can  be  repeated  for  both  the  x-  and  y- tilt  channels  for  all  69  illu¬ 
minated  subapertures,  continuing  with  the  noise-free  measurement  assumption.  The  138 
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resulting  equations  can  be  augmented  into  a  matrix  form: 


£(*>)=  ^Na(!;)  (3.49) 

Therefore,  the  H'  matrix  of  Equation  (3.39)  is: 

H'  =  (3.50) 

The  Fortran  program  developed  to  calculate  the  N  matrix  for  this  research  is  archived 
at  AFIT  [39].  Appendix  F  shows  the  resulting  N  matrix. 

To  determine  the  value  of  the  scaling  constant,  L/A,  results  of  previous  AFIT  research 
are  useful.  Miller  [27]  collected  measurement  data  from  a  realization  of  the  Hartmann 
sensor  studied  in  this  research.  He  excited  the  sensor  with  543.5-nm  laser  light  and  recorded 
the  outputs  from  the  slope  channels.  These  outputs  were  not  in  wavelengths  of  tilt,  but 
rather  in  internal  programmable  microcoded  processor  (PMP)  units.  He  determined  that 
one  wavelength  of  tilt  across  a  subaperture,  on  the  average,  resulted  in  a  slope  channel 
output  of  745  PMP  units.  Simple  algebra  shows  that  one  wavelength  of  y-tilt  across  a 
subaperture  corresponds  to  the  Zernike  coefficient  a-i  —  2.574.  The  value  of  L/A  applicable 
to  his  research  can  be  determined  by  considering  the  non-zero  elements  of  the  aj  column 
of  the  N  matrix,  which  are  2.0419.  Solving  the  equation: 

745  =  (2.0419)  (2.574)  (3.51) 

yields  a  value  for  of  141.75. 

To  keep  the  analysis  as  general  as  possible,  this  research  assumes  the  outputs  of  the 
Hartmann  sensor  are  the  x-  and  y-tilts  of  each  subaperture,  in  units  of  514-nm  wavelengths 
( not  PMP  units).  Similar  algebra  yields  a  value  for  j  of  0.19025.  The  N  matrix  derived 
in  Appendix  F  should  be  premultiplied  by  this  constant  to  get  the  H’  matrix  used  in 
Equations  (3.39). 
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5\S.o  Derivation  of  R  T  he  outputs  of  wavefront  sensors  are  generally  corrupted  by 
noise.  Numerous  literature  references  [30:139]  [14:30]  [40:1772]  [40:1773]  [32:821]  indicate 
the  noise  is  zero-mean,  white,  Gaussian,  unccrrelated  from  subaperture  to  subaperture, 
and  uncorrelatea  o  the  wavefront  phase.  Possible  sources  of  this  noise  include  photon 
shot  noise,  A/D  quantization  noise  (assuming  the  Hartmann  outputs  are  digital),  and 
thermal  noise.  The  photon  shot  noise  variance  is  inversely  proportional  to  the  number  of 
photons  counted  during  a  measurement  cycle  [47],  The  photon  count  is,  in  turn,  related  to 
light  intensity.  The  lower  the  light  intensity,  the  more  shot  noise  there  will  be  in  the  slope 
measurement.  A/D  quantization  noise  is  typically  uniformly  distributed  between  minus 
1/2  and  plus  1/2  least  significant  bit.  This  distribution  could  be  approximated  as  Gaussian 
with  the  3-c  value  set  to  1/2  the  least  significant  bit  [24:364].  Thermal  noise  is  related  ti> 
absolute  temperature,  through  Boltzman’s  constant. 

Since  the  intensity  of  the  viewed  object  is  expected  to  vary  greatly  [23]  as  the  target 
changes  orientation  with  respect  to  the  viewer,  and  also  from  target  to  target,  it  is  likely 
that  the  shot  noise  contribution  to  the  measurement  noise  will  also  vary.  This  research 
treats  such  noise  as  the  dominant  noise  source,  ignoring  all  others  [47],  Therefore,  this 
’•csearch  includes  a  parameter  study  on  the  R  matrix.  Using  equations  discussed  in  Welsh 
and  Gardner  [47],  this  research  calculates  slope  noise  strengths  corresponding  to  1000,  100, 
and  10  photons  per  subaperture.  Table  3.2  shows  values  used  for  the  diagonal  elements  of 
the  R  matrix. 


Table  3.2.  Measurement  Noise  Strengths 


N 

DESIGNATOR 

Ra  (wavelength2) 

1000 

Low 

0.01755 

100 

Medium 

0.0555 

10 

High 

0.1755 
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S.G  Summary 

This  research  develops  stochasti-'  models  for  atmospheric  turbulence  <  Tects,  de¬ 
formable  mirror  dynamics,  and  a  Hartmann-type  wavefront  sensor.  The  phase  distortion 
due  to  atmospheric  turbulence  is  modeled  in  the  Zernike  functional  space.  Tne  Zernike 
coefficients  are  modeled  as  first-order  Gauss-Maikov  random  processes.  The  models  are 
obtained  by  fitting  theoretically-derived  autocorrelation  data  to  appropriate  functions. 
Using  data  from  the  AFWL,  this  research  synthesizes  a  deterministic  dynamics  model  for 
a  deformable  mirror.  The  ability  of  the  mirror  to  reconstruct  the  14th  Zernike  function 
is  verified  in  a  simple  simulation.  The  measurement  matrix  H  is  calculated  by  projecting 
each  Zernike  mode  of  distortion  into  subaperture  slope  space,  incorporating  an  appropri¬ 
ate  weighting  function.  Finally,  values  of  slope  sensor  noise  due  to  phoion  shot  noise  are 
calculated,  to  be  used  as  a  varying  parameter  in  simulations. 
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IV.  Controller  Design 


4-1  Kalman  Filter  Design 

The  Kalman  filter  for  this  research  is  ef  the  same  order  as  the  truth  model.  Since  both 
the  atmosphere  and  the  mirror  Zexnike  coefficients  are  first-order  processes,  the  Kalman 
filter  is  a  28- state  filler — 14  for  the  atmosphere  and  14  for  the.  mirror.  For  the  time-invariant 
model  assumed  here,  the  continuous-time  dynamics  and  discrete-time  measurement  truth 
model  equations  are: 


x(0 

=  Fx(0  +  Bu(f)  +  w/0 

(4.1) 

&(*;) 

=  Hx(ti)Ty(t,) 

(4.2) 

Tne  r  anu  Q  matrices  are  diamond,  With  elements  given  m  Tabic  -1.1.  Tire  y  tilt 
and  x-tilt  no:se  strengths  have  bee  r  attenuated  to  reflect,  the  presence  of  tilt  mirrors. 
These  tilt  mirrom  are  modeled  as  removing  95  percent  of  atmosphere-induced  tilt.  The 
B  matrix  is  merely  determined  using  Equations  (3.35)  and  (3.37);  the  M  matrix  is  given 
i'i  Appendix  E.  The  H  matrix  is  determined  from  Equations  (3.39),  (3.50),  and  the  N 
matrix  of  Appendix  F.  The  noise  variance  matrix,  R  is  an  identity  matrix  premultiplied 
by  one  of  the  values  of  Table  3.2,  depending  on  the  noise  model. 

The  filter  models  the-  continuous-time  dynamics  of  the  truth  model  in  the  discrete 
domain  as  suggested  by  Equation  (2.38),  forgoing  on-line  integration.  Table  4.2  summarizes 
the  dimensionality  of  the  discrete  state-space  filter  model  for  the  adaptive  optics  system. 

With  a  sample  period  of  7  milliseconds  (the  maximum  sampling  rate  of  the  Reticon 
chip  in  the  Hartmann  sensor  [27]),  the  state  transition  matrix  for  this  time-invariant  system 
model  is  calculated  as  [24:42]: 


-  fl)  =  e(F)(°.007) 
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(4.3) 


Table  4.1.  Continuous  State-Space  Model 


(wavelenglh2/second2) 

0.1605 

0.075 

0.405 

0.800 

0.509 

0.124 

0.223 

0.210 

0.241 

n  o/'o 
O.oui) 

0.0694 

0.0833 

0.0645 

0.117 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 
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Table  4.2.  Dimensionality  of  Kalman  Filter  Vectors  and  Matrices 


Matrix/Vector 

Dimension 

28  X  28 

Z.(U) 

28  X  1 

Bd 

28  X  97 

Mb) 

97  x  1 

2£(*i) 

28  X  1 

Qd 

28  x  28 

z(t.) 

138  X  1 

H 

138  x  28 

v(b) 

138  x  1 

R 

138  x  138 

Since  F  is  a  diagonal  matrix,  so  is  <F(t,+1  -  U),  and  the  exponentiation  is  easily 
accomplished  term-by-term.  The  discrete  version  of  the  dynamics  driving  noise  is  deter¬ 
mined  using  Equation  (2.39).  For  this  time-invariant  model,  with  G  being  the  identity 
matrix,  and  F  being  diagonal,  Equation  (2.39)  is  manipulated: 


Qd  = 


i  +  l  n n  'T 

/  0(t,+i  -  r)GQGi<l>1(t,+1  -  r)  dr 

Jt, 

r  o.oor  _ 

_  /  eF(0.007-rjQeF(0.007-T) 

Jo 

/•O.U07 

=  Q  ^(0.007-r)  dr 

JO 


_  f  0.007 

=  Qt  2F0.007  /  e-2Fr 

Jo 


=  Qe2F0007  L-2Ft 


-F 


-l 


.007 


_  q62F0.0O7  |e-2F0.007  _  jj  |  Tlp-1 


Q  7—1  J  _  e-2F0.007 


(4-4) 
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The  diagonal  elements  of  —  t,)  and  Qj  are  given  in  Table  4.3.  Non-diagoi:  ..  elements 

are  zero. 

The  Kalman  filter  processes  the  138  slope  measurements  every  0.007  seconds  using 
Equations  (2.43)  through  (2.47).  The  filter  assumes  the  measurements  are  available  in¬ 
stantly.  The  important  output  of  the  filter  is  the  estimate  of  the  system  slate,  since  this 
is  ultimately  multiplied  by  the  controller  gain. 

The  Kalman  filter  requires  initial  conditions  on  both  the  state  estimate  x(0)  and 
covariance  P(0).  Since  the  states  model  the  Zernike  coefficients  as  zero-mean,  the  initial 
filter  estimate  of  the  system  state  is  selected  to  be  a  zero  vector.  The  first  14  diagonal 
elements  of  P(0)  are  set  to  large  values  to  rellect  initio!  uncertainty. 

4-2  Linear  Quadratic  Regulator  Design 

Equation  (2.41)  defines  the  cost  to  be  minimized  by  the  quadratic  regulator.  For  the 
adaptive  optics  system  of  interest,  it  makes  sense  to  assign  cost  to  the  phase  distortion  in 
the  reflected  image,  since  this  is  the  image  which  is  desired  to  be  distortionless.  It  also 
makes  sense  to  assign  cost  to  the  control  voltages,  since  they  are  restricted  to  be  within  a 
finite  range  (±10  volts).  The  derivation  of  both  the  X  and  the  U  cost  matrices  are  now 
discussed. 

The  phase  distortion  in  the  reflected  image  is  modeled  by  its  set  of  Zernike  coefficients. 
This  set  of  Zernike  coefficients  is  obtainable  from  the  discrete  state  vector  by: 


The  quadratic  cost  associated  with  the  distoition  in  the  reflected  image  is: 
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Table  4.3.  Filter  Stare  Transition  Matrix  and  Discrete  Driving  N 


j 

^1  +  1  •)  t  t  ) 

Qdjj  (wavelength2) 

1 

0.904810273562201 

0.004073001858519 

2 

0.954383210022565 

0.002004933000163 

3 

0.920977197358095 

0.010455681053018 

4 

0.737595286709826 

0.0 1C7784 06543654 

5 

0.818738942175464 

' 

0.011746606438655 

. 

wm 

0.823279104805153 

0.00287G4 74222 156 

n 

0.704688089718714 

0.004490459090181 

8 

0. 737595286 709S26 

0.004404331717709 

9 

0.7046SS08971S714 

0.004852917671451 

10 

0.943357461794596 

0.0095937 196 16645 

11 

0.677788484408774 

0.001350533885479 

12 

0.627118334037179 

0.001516351317027 

13 

0.900765787779704 

0.001629746010087 

14 

0.627074461402559 

0.002129683020039 

15 

0.000000175785735 

0.000000000000000 

1G 

0.000O00175785735 

0.000000000000000 

17 

0.0000001 75785735 

0.000000000000000 

IS 

0.000000175785735 

0.000000000000000 

SB 

0.000000175785735 

0.000000000000000 

m 

0.000000175785735 

0.000000000000000 

21 

0.000000175785735 

0.000009000000000 

22 

0.000000175785735 

0.000000000000000 

23 

0.000000175785755 

0.000000000000000 

24 

0.000000175785735 

0.000000000000000 

25 

0.000000175785735 

0.000000000000000 

20 

0.000000175785735 

0.000000000000000 

27 

0.000000175785735 

0.000000000000000 

2S 

0.000000175785735 

0.000000000000000 
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where  the  C  is  a  diagonal  matrix  of  cost  elements.  Assume  one  desires  Hie  maximum  rms 
phase  deviation  due  to  any  Zernike  mode  to  be  1/20  of  a  wavelength.  Remembering  that, 
the  Zernike  coefficients  map  directly  into  rms  phase  distortion,  one  can  say  the  maximum 
value  for  any  Zernike  coefficient  is  1/20.  Each  diagonal  term  of  C,  for  the  first  design,  can 
be  set  equal  to  the  reciprocal  of  the  square  of  the  maximum  value  [20:69],  Therefore,  the 
diagonal  dements  of  C  are  400,  with  all  other  elements  zero,  for  initial  tuning. 

Manipulating  the  expression  in  Equation  (4.6)  one  obtains: 


xJ(U)  x£(*i) 


ATCA  atc 
CA  C 


Zoiu) 

25m  (<i) 


which  is  exactly  of  the  form 

x^X  x 


(4.7) 


(4.8) 


The  extraction  matrix  A  is  a  14-by-14  identity  matrix  for  this  research,  since  each  atmo¬ 
spheric  Zernike  coefficient  is  modeled  as  first-order  Gauss- Markov.  Therefore,  the  state 
cost  matrix  X  is  a  28-by-28,  singular,  banded  diagonal  matrix: 


X  = 


c  c 
c  c 


(4.9) 


The  first  design  iteration  of  the  cost  matrix  U  associated  with  the  control  voltages 
is  simpler  to  derive.  The  maximum  allowed  control  voltage  to  any  actuator  is  ±  10  volts. 
Using  the  reciprocal- squared  method  results  in  the  initial  design  for  the  U  matrix  being 
a  diagonal  matrix  with  diagonal  elements  of  0.01.  If  the  simulated  peiformance  of  the 
adaptive  optics  system  is  unacceptable,  these  cost  matrices  can  be  varied  to  try  to  improve 
performance. 

Knowing  the  <£>((, +1  —  Be,  X,  and  U  matrices,  one  can  use  the  Martix-X  command 
“DREGULATOR”  to  obtain  the  steady-state  solution  of  the  backward  Riccati  equation 
and  obtain  the  optimal  steady-state  controller  gain  matrix  G*  (See  Equations  (2.51)  and 
(2.52)).  This  value  of  G'  can  now  be  used  as  the  controller  gain  in  a  digital  simulation 
of  the  adaptive  optics  system.  The  cost  matrices  X  and  U  can  be  adjust*  ,!  to  try  to 
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improve  performance.  The  controller  multiplies  its  optimal  gain  by  the  state  estimate 
from  the  Kalman  filter  to  obtain  a  set  of  control  voltages  for  the  mirroi  actuators  (See 
Equation  (2.50)). 

4.3  Summary 

This  chapter  assembles  the  results  of  previous  chapters  into  suitable  Kalman  filter 
and  controller  form.  The  cost  matrices  X  and  U  are  calculated  for  the  first  iteration  of 
controller  design.  Matrix-X  software  is  used  to  solve  for  the  steady-state  controller  gain. 
Depending  on  simulation  results,  the  cost  matrices  may  need  to  be  tuned  to  maximize 
performance,  i.e.,  minimize  rms  phase  distortion  in  the  reflected  image. 
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I  .  Svmilot  ton  Jicsnltf 


5.1  Methodology 

A  digital  simulation  of  the  adaptive  optics  system  was  desired  in  order  to  test  the 
effectiveness  of  the  LQG  controller  in  reducing  phase  aberration.  The  Multimode  Simula¬ 
tion  for  Optimal  filter  Evaluation  (MSOFE)  software  tool  was  used  [28],  The  simulation 
involved  implementing  a  truth  model  of  the  atmosphere/mirror/sensor  system  to  generate 
simulated  discrete-time  Hartmann  sensor  measurements.  The  Kalman  filter,  also  em¬ 
bedded  in  the  simulation,  processed  these  simulated  measurements  and  determined  the 
estimate  of  the  28  truth  model  states.  The  LQ  regulator  then  multiplied  this  state  vector 
estimate  by  the  optimal  controller  gain  (pro-calculated  ofT-line),  and  the  resulting  controi 
voltage  commands  were  sent  to  the.  simulated  mirror  actuators.  The  nonlinearity  of  the 
mirror  response  to  applied  control  voltage  was  simulated  by  limiting  the  mirror  voltage  to 
the  ±  10  volt  region.  For  example,  if  the  controller  tried  to  command  an  11 -volt  control 
voltage  to  an  actuator,  logic  in  the  simulation  program  truncated  the  value  to  10  volts. 
This  nonlinearity  imposed  the  use  of  multiple  Monte  Carlo  simulations,  as  opposed  to  a 
single  covariance  analysis  [24:329]  [2S]. 

Outpu  s  from  the  simulation  included  time  histories  of  the  28  truth-model  states, 
the  filter's  estimates  oi  these  states,  the  filter  covariance  matrices,  and  the  maximum  and 
minimum  control  voltages  of  the  mirror  actuators.  Post-simulation  data  reduction  yielded 
time  histories  of  the  nns  phase  distortion  at  the  optics  entrance  aperture  as  well  as  the  rms 
phase  distortion  after  correction.  The  Matrix- X  software  package  was  used  to  perform  the 
data  reduction  as  well  as  generate  plots  of  time  histories.  In  addition  to  single  realizations 
of  these  time  histories,  ensemble  statistics  were  generated  using  Monte  Carlo  analysis  of 
ten  runs. 

As  previously  discussed,  there  is  expected  to  be  significant  variation  of  measurement 
noise  strength,  i.e.,  the  elements  of  the  R  matrix,  in  a  real  adaptive  optics  system.  The 
simulation  was  performed  in  nine  “studies”  in  order  to  investigate  performance  sensitivity 
to  R.  Each  study  simulated  a  difh  rent  combination  of  truth  model  and  filler  model 
measurement  noise  strengths.  Table  5.1  shows  the  combinations  simulated. 
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Table  5.1.  Truth  and  Filter  Model  Measurement  Noise  Strengths 


STUDY 

TRUTH 

FILTER 

1 

Low 

Low 

2 

Medium 

Medium 

3 

High 

High 

4 

Low 

Med 

5 

Low 

High 

C 

Medium 

Low 

7 

Medium 

High 

8 

High 

Low 

9 

High 

Medium 

The  noise  strengths  of  “Low”,  “Medium”,  and  “High”  correspond  to  photon  counts 
per  subaperture  of  1000,  100,  and  10,  respectively;  see  Table  3.2.  The  objective  here  was 
to  determine  the  effect  of  mismodeling  the  noise  strength  in  the  Kalman  filter.  In  a  real 
implementation,  R  in  the  filter  model  may  be  fixed  and  therefore  wrong  for  the  case  of 
the  true,  intensity-dependent  measurement  noise.  Since  the  regulator  design  was  based  on 
a  deterministic  version  of  the  stochastic  dynamic  equations  due  to  certainty  equivalence, 
changing  the  R  matrix  does  not  affect  the  steady-state  value  of  the  regulator  gain  matrix 
G*.  Thus,  all  nine  studies  used  the  same  regulator  gain  matrix.  The  following  discussion 
presents  a  verbal  and  graphical  description  of  study  1. 

5.2  Study  J  Description 

Study  1,  having  low  measurement  noise  cor  ?c.tly  modeled  in  the  Kalman  filter,  is 
the  most  optimistic  of  the  studies  from  a  performance  perspective.  Since  the  mirror  (states 
15-2.8)  is  assumed  deterministic  in  both  the  truth  and  filter  models,  the  filter  estimates  of 
the  mirror  states  are  trivial,  and  they  will  be  omitted  from  most  of  the  discussion.  In  an 
operational  system  these  deterministic  states  should  rot  even  be  included  in  the  Kalman 
filter,  for  computational  reasons.  Their  inclusion  in  the  filter  here  is  merely  convenient. 


STATE  XI  (WAVES) 


5.2.1  Simulated  Atmospheric  State  Behavior  Th;  research  models  image  phase  dis¬ 
tortion  induced  by  atmospheric  turbulence  as  fourteen  time- varying  Zernike  coefficients. 
It  is  arbitrarily  assumed  that  tilt  mirrors  ahead  of  the  adaptive  optics  system  remove  9 5 
percent  of  the  first  two  Zernike  modes  (y  -  tilt  and  x-tilt).  Zernik''  coefficients  for  the  re¬ 
maining  tilt  distortion  as  well  as  the  other  twelve  modes  are  states  1-14  of  the  truth  model. 
Figure  5.1  shows  a  time-history  of  a  sample  realization  of  the  y-tiit  Zernike  coefficient  x\ 
as  well  as  the  filter  estimate,  x;. 


Figure  5.1.  Atmospheric  Y-tilt  State  and  Filter  Estimate 

Filter  error  can  be  defined  as  the  true  state  minus  the  filler  estimate  of  the  state. 
The  filter  covariance  P  is  the  filter’s  indication  of  uncertainty  in  its  estimates,  as  in  liqua¬ 
tion  (2.42).  The  square  root  of  the  (1,1)  element  of  the  P  matrix  is  what  the  filter  believes 
to  be  the  ]-o  value  of  its  error.  Although  Figure  5.1  shows  that  the  filter  appears  to  be 
tracking  the  true  state  value,  it  docs  not  indicate  how  the  actual  filter  error  compares  with 


what  it  thinks  are  its  1  -a  values.  Figure  5.2  shows  a  time  history  of  the  filter  error  for 
state  1,  as  well  as  these  filter-computed  1  -a  values. 


0  .05  .1  .15  .2  .25  .3  .35  .4  .45  .5 
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Figure  5.2.  Atmospheric  Y-tilt  Filter  Error  and  Filter  Variance 


Looking  at  the  single-sample  realization  of  filter  error  of  Figure  5.2  and  assuming  the 


error  is  an  ergodic  random  process,  it  appears  that  1]  is  a  reasonable  l-</  value 

for  the  error,  and  that  the  error  is  zero  mean.  In  order  to  obtain  the  “trac”  error  mean 
and  variance,  a  Monte  Carlo  analysis  of  ten  runs  was  accomplished.  The  equations  used 
to  process  the  results  of  these  ten  runs  are: 


1  50 

a*e(0  ~ 


10 


LTE(*l«>-».fc«>j 
*=  2 


(5.1) 

(5.2) 


C£(0  ~ 


where 


x*;(f)  -  Xfc(<)  (wavelengths) 
sample  realization  number 
mean  of  random  process  E(<)  (wavelengths) 
variance  of  random  process  E(t)  (wavelengths') 

Figure  5.3  shows  a  plot  cf  the  mean  and  standard  deviation  of  the  filter  error  for 
state  1,  as  calculated  from  the  ten  sample  realizations  using  Equations  (5.1)  and  (5.2). 
Visual  inspection  of  this  plot  reveals  the  error  process  is  approximately  zero  mean,  with 
standard  deviations  approximately  equal  to  the  yOhl,  l)  values  from  Figure  5.2.  If  more 
sample  realizations  had  been  included  in  the  Monte  Carlo  analysis,  the  similarity  would 
most  likely  be  even  greater. 

Plots  similar  to  those  of  Figures  5.1,  5.2.  and  5.3,  can  be  found  in  Appendix  G  for 
the  remaining  13  states  corresponding  to  image  phase  distortion  prior  to  correction.  This 
complete  set  of  plots  is  for  study  1  only. 

b.2. 2  Performance  An.  lysis  The  performance  of  the  LQG  ontroller  can  be  ex¬ 
pressed  in  terms  of  rms  phase  distortion  in  the  corrected  image  versus  the  ::ms  phase 
distortion  of  the  incident  image.  Again,  incident  here  means  after  the  tilt  mirrors  have 
removed  most  of  the  gross  tilt..  Slates  1 — 14  of  the  system  model  are  the  Zernike  coeffi¬ 
cients  for  atmosphere-induced  phase  distortion,  and  states  15 — 28  are  Zernike  coefficients 
for  the  minor-induced  “counterdislortion” .  Since  reflection  from  tbs  deformable  mirror  is 
modeled  as  addition  of  atmospheric  and  minor  phase  distortions,  the  Zernike  coefficients 
cf  the  reflected  (i.e.,  corrected)  image  are: 
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Figure  5.3.  Mean  and  Standard  Deviation  of  Atmospheric  Y-tilt  Filter  Error 


It  has  already  been  shown  that  the  rms  phase  distortion  is  the  square  root  of  the 
sum  of  the  squares  of  the  Zernike  coefficients  (see  Appendix  C).  The  rms  phase  error  of 
the  incident  image  is: 


14 


and  for  the  corrected  image,  the  rms  phase  distortion  is: 


(5.6) 


C^rmj(0 


14 


N 


X>,  +  x!+14): 

>=i 


(5.7) 


A  time  history  of  a  sample  realization  of  the  rms  phase  distortion  is  shown  in  Figure  5.4 
for  both  the  incident  and  corrected  images.  The  values  for  this  plot  were  obtained  during 
the  post-simulation  data  reduction  phase  of  the  study,  using  Equations  (5.3)  and  (5.7). 
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T lie  plot  clearly  shows  the  reduction  in  rms  phase  error.  Moreover,  this  single  realization 
shows  the  rms  phase  error  of  the  corrected  image  tends  to  stay  in  tln>  vicinity  of  nix  ■  6.1 
wavelen  :,hs,  whereas  the  rms  phase  error  or  li-c  incident  imago  varies  between  about  0.1 
wavelength  and  0.5  wavelength. 

Calculating  similar  results  for  each  of  t  he  ten  Monte  Cario  realizations  and  gene;  ating 
the  usual  statistics  results  in  values  plotted  in  Figure  5.5.  The  upper  three  lines  of  the 
graph  show  the  mean  and  the  [meanil-cr]  values  on  the  rms  phase  distortion  of  the  incident 
image.  The  apparent  transient  in  these  three  lines  during  the  initial  0.05  seconds  of  the 
simulation  was  caused  by  unrealistic  initial  condition'  in  the  truth  model  (.all  dates  zero). 
It  would  be  a  more  realistic  simulation  if  the  initial  true  states  were  random.  The  lower 
three  lines  of  the  graph  represent  the  mean  and  the  [mean  ±  \-c)  val  ms  on  tlm  uns  phase 
distortion  of  the  corrected  image.  One  feature-  of  note  is  that  the  l-e  values  on  the  corrected 
image  are  tighter  than  on  the  incident  image.  This  seems  to  suggest  the  quality  of  the 
corrected  image  is  somewhat  constant,  despite  wide  variation  in  the  amount  of  atmospheric 
distortion.  Another  feature  is  the  sa.wlooth  appearance  of  the  lower  set  of  plots.  The  period 
of  the  sawtooth  appears  to  be  the  Hartmann  sampling  period  fO.OOT  sec).  This  suggests 
the  correction  is  most  effective  just  after  a  measurement,  and  degrades  as  the  atmosphere 
changes  between  measurements.  One  may  wonder  if  even  better  performance  m  possible. 
In  other  words,  one  can  ask,  what  is  limiting  the  performance  shown  in  Figure  5.5?  At 
least  three  answers  are  possible:  1)  saturation  of  mirror  actuators,  2)  improper,}'  chosen 
weighting  matrices  in  the  LQ  regulator,  or  3)  the  filler's  estimation  errors.  The  first 
possible  reason,  actuator  saturation,  is  immediately  ei'iumated  based  on  Figure  5.6.  This 
plot  shows  the  absolute  envelope  for  the  actuator  control  voltages  for  tiie  ensemble  of 
ten  Monte  Carlo  realizations.  Given  that  saturation  occurs  when  the  magnitude  of  che 
commanded  voltage  exceeds  10  volts,  and  the  maximum  actual  excursion  was  only  about 
±  one  volt,  saturation  did  not  occur.  In  fact,  one  could  say  the  mirror  had  quite  a  bit 
of  ‘'remaining”  capability  left.  The  second  reason,  improperly  chosen  weighting  matrices 
X  and  U,  was  eliminated  b}  a  tuning  experiment.  The  nonzero  elements  of  the  X 
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Figure  5.6.  Control  Voltage  Envelope  for  Monte  Carlo  Study 


matrix  of  the  cost  equation  were  increased  by  50  percent.  This  penalized  the  distortion  of 
the  corrected  image  more  heavily.  The  simulation  was  re-run  using  the  new  steady-state 
controller  gain,  G’,  and  performance  did  not  improve.  The  third  reason,  filter  estimation 
error,  was  analyzed  in  the  following  manner.  The  filter’s  function  is  essentially  to  estimate 
the  Zernike  coefficients  of  the  incident,  uncorrected  image.  Assume  the  regulator/mirror 
combination  can  perfectly  implement  the  filter’s  estimate  (without  a  sign  change).  The 
rms  phase  error  due  solely  to  the  filter’s  estimation  error  can  be  thought  of  as  the  lower 
bound  on  rms  phase  error  attainable  and  is  calculated  using: 


vVm  j  ( l  ) 


14 


^Dx>  x>)? 


(5.3) 


Monte  Carlo  analysis  of  values  calculated  using  Equation  (5.3)  is  shown  in  Figure  5.7.  The 
extreme  similarity  of  this  graph  and  the*lower  three  plots  of  Figure  5.5  strongly  suggest 
that  filter  estimation  error  is  the  performance-limiting  factor.  An  alternate  method  of 
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determining  the  contribution  of  filter  estimation  error  to  the  rms  phase  distortion  would 
be  to  input  x  instead  of  x  into  the  controller,  and  see  if  the  rms  error  is  eliminated.  This 
alternate  method  was  not  used  in  this  research. 


—  Mean 
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Figure  5.7.  Monte  Carlo  Study  of  RMS  Phase  Distortion  Caused  by  Filter  Error 


5.3  Performance  Results 

The  above  analysis  section  dealt  exclusively  with  study  1,  the  case  of  low  measure¬ 
ment  noise  correctly  modeled  in  the  filter.  Similar  graphs  for  the  remaining  eight  Studies 
are  shown  in  Appendix  K.  ,'Also  shown  are-  filter  estimation  plots  for  Zernike  modes  1 
and  14.)  In  order  to  compare  results  from  the  nine  studies  conveniently,  further  data 
compres.  1  was  accomplished.  This  amounted  to  time-averaging  the  mean-value  plots  of 
rms  phase  distortion  from  each  of  the  nine  Monte  Carlo  studies.  Thus  each  of  the  nine 
studies  is  summarized  by  three  numbers:  the  rms  phase  distortion  of  the  incident  image, 
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the  rms  phase  distortion  of  the  corrected  image,  and  the  rms  phase  distortion  which  could 
be  achieved  if  the  regulator/mirror  could  perfectly  implement  the  filter’s  estimate,  i.e.,  the 
filter  “error.”  Furthermore,  all  three  of  these  numbers  can  be  normalized  by  dividing  by 
the  first.  This  normalization  provides  a  fairer  comparison  across  the  nine  studies.  Table  5.2 
shows  the  results  of  this  time- averaging  and  normalization. 


Table  5.2.  Summary  of  Adaptive  Optics  System  Simulated  Performance 


STUDY 

RMS 

Inc. 

DISTORT 

Cor. 

ION 

Fit. 

Inc. 

s'ORMALlZ 

Ccr. 

ED 

Fit. 

1 

0.2552 

0.1083 

0.1031 

1 

0.4244 

0.4040 

2 

0.2557 

0.1126 

0.1079 

1 

0.4404 

0.4220 

3 

0.2598 

0.1321 

0.1285 

1 

0.5085 

0.4946 

4 

0.2546 

0. 1088 

0.1033 

* 

1 

0.427 J 

]£^|||j9[ 

5 

0.2536 

0.1112 

1 

0.46S1 

0.4385 

6 
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As  expected,  study  1  resulted  in  the  lowest  rms  phase  distortion  in  the  corrected 
image.  This  was  the  case  oflow  measurement  noise,  correctly  modeled  in  the  filter.  That 
the  rms  phase  distortion  of  the  corrected  image  was  only  slightly  larger  than  the  rms  phase 
error  caused  by  the  filter  error  indicates  that  most  of  the  phase  error  in  the  corrected  image 
was  due  to  filter  estimation  errors.  As  the  true  measurement  noise  increased  (studies  1  — > 
2  —  3),  and  the  filter  R  was  modified  correspondingly,  the  performance  was  progressively 
poorer.  This  reflects  the  fact  that  noisier  measurements  resulted  in  a  less  accurate  state 
estimation,  even  with  the  filter  properly  tuned. 

The  worst  performance  came  from  sV'dy  8.  This  was  the  case  of  high  measurement 
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noise,  mismodeled  in  the  filter  as  low  measurement  noise,  lleuristically,  the  filter  had  too 
much  confidence  in  the  incoming  measurements,  and  tended  to  disregard  it’s  own  model  of 
the  dynamics.  The  result  was  a  corrected  image  having  an  rms  phase  distortion  of  0.1428 
wavelengths. 

5.4  Summary 

This  chapter  discusses  the  digital  simulation  of  the  adaptive  optics  system.  The 
MSOFE  software  [28]  is  used  to  accomplish  the  simulation.  This  software  simulates  both 
real-world  behavior  and  Kalman  filter  processing.  Modifications  to  the  software  allow 
for  implementation  of  LQG  control.  The  simulations  comprise  a  set  of  nine  studies,  each 
having  a  different  combination  of  true  measurement  noise  and  filter  model  thereof.  Study 
1,  the  case  of  low  measurement  noise  correctly  modeled  in  the  filter,  is  analyzed  in  detail. 
A  table  of  results  is  presented  which  indicates  filter  estimation  errors  limit  the  controller  s 

performance. 
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VI.  Conclusions  and  Recoin  men  dal  ions 

6.1  Summary 

This  research  considered  the  design  of  a  nominal  linear  quadratic  Gaussian  (LQG) 
controller  for  a  ground-based  adaptive-optics  telescope.  Phase  distortion  caused  by  at¬ 
mospheric  turbulence  was  modeled  as  14  time-varying  Zcrnike  coefficients.  The  effect  of 
tilt  mirrors  was  modeled  as  removing  95  percent  of  the  mean  square  contribution  of  the 
first  two  Zcrnike  modes.  Dynamics  of  the  97-acluator  deformable  mirror,  following  the  lilt 
mirrors  in  the  optical  path,  were  modeled  as  a  14  deterministic  first-order  lags. 

A  69-subuperture  Hartmann-type  wavefront  sensor  was  assumed  to  be  the  measure¬ 
ment  device.  The  slope  outputs  from  the  subapertures  comprised  a  138-element  measure¬ 
ment  vector.  The  sampling  period  was  assumed  to  be  7  milliseconds,  corresponding  to 
the  maximum  detector  rate.  A  2S-stale  Kalman  filter  processed  the  measurenm.  ts  from 
the  wavefront  sensor  and  obtained  estimates  of  system  states.  A  constant-gain  linear 
quadratic  (LQ)  regulator  processed  these  stale  estimates  and  determined  an  appropriate 
set  of  commands  for  the  deformable  mirror. 

The  entire  control  system  was  digitally  simulated  using  the  Multimode  Simulation 
for  Optimal  Filter  Evaluation  (MSOFE)  software.  Nine  simulation  studies  were  conducted 
to  investigate  the  effects  of  mismodeling  the  noise  in  the  measurement  device. 

6.2  Conclusions 

1.  The  LQG  approach  used  in  this  research  makes  sense  because  the  desired  goal  of 
reducing  the.  rms  phase  error  translates  directly  into  the  quadratic  cost  criteria.  The 
mean-square  phase  distortion  is  the  sum  of  the  squares  of  the  Zcrnike  coefficients 
(states). 

2.  Most  of  the  open  literature  models  atmospheric  turbulence  as  having  Kolmogorov 
properties.  Taylor’s  frozen  field  assumption  is  also  popular.  Results  based  on  actual 
measurements  of  atmospheric  effects  on  image  quality  are  sparse. 
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3.  Assuming  the  atmospheric  distortion  is  well-modeled  by  the  Gauss-Murkov  processes 
of  this  research,  the  deformable  mirror  operates  in  the  linear  region,  i.e.,  its  actuators 
are  seeing  control  voltages  in  the  ±1  volt  range,  whereas  saturation  occurs  at  the 
±10  volt  limits. 

4.  Modifications  to  the  MSOFE  software  allowed  for  simulation  of  LQG  control  These 
consisted  of  pre-multiplication  of  the  filler  state  estimate  vector  (after  measurement 
update)  by  the  steady-state  controller  gain,  pre-multiplication  of  this  result  by  the  in¬ 
put  distribution  matrix,  and  adding  the  result  to  the  right-hand  side  of  the  dynamics 
equation. 

5.  Based  on  Monte  Carlo  simulation  results,  the  adaptive  optics  system  did  reduce 
phase  distortion  in  all  cases.  Table  5.2  shows  that,  the  (simulated)  deformable  mirror 
reduces  the  rms  phase  distortion  to  about  40 — GO  percent  of  its  incident  (post-tilt 
mirror)  value.  The  best  performance  was  achieved  when  the  filter  model  of  the 
measurement  noise  was  “correct”  relative  to  the  truth  model.  When  the  measurement 
noise  in  the  truth  model  increased,  the  performance  degraded,  even  if  the  filter  was 
retuned  correspondingly. 

0.  A  single  realization  of  the  performance  history  (Figure  5.4)  hints  that  the  phase  qual¬ 
ity  of  the  corrected  image  is  somewhat  constant,  regardless  of  the  actual  magnitude 
of  the  incoming  atmospheric  phase  distortion.  This  result  is  most  likely  related  to 
the  operation  of  the  actuators  far  from  their  saturation  limits. 

T.  Comparison  of  the  remaining  phase  distortion  after  correction  to  the  phase  distortion 
caused  by  filter  estimation  error  (Figures  -.5  and  5.7)  indicates  that  state  estimation 
error  is  the  factor  most  limiting  performance. 

C.  5  llccrmvucndatwnt 

C.’J.l  Modeling  This  researcli  is  based  on  a  set  of  nominal  models,  mostly  derived 
front  theoretical  results.  The  aforementioned  “performance’'  of  the  LQG  control  law  is 
only  as  valid  as  the  truth  model  of  the  real  world.  This  research  never  claimed  to  develop 
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